#Loading required libraries

library(phyloseq)
# Define packages
## CRAN repositories
cran_packages <- c("bookdown", "knitr", "tidyverse", "plyr", "grid", "gridExtra", "kableExtra", "xtable", "ggpubr")
## Bioconductor repository
bioc_packages <- c("phyloseq", "dada2", "DECIPHER", "phangorn", "ggpubr", "BiocManager","DESeq2", "microbiome", "philr")
##  GitHub repository
git_source <- c("twbattaglia/btools", "gmteunisse/Fantaxtic", "MadsAlbertsen/ampvis2", "opisthokonta/tsnemicrobiota") # fuente/nombre del paquete
git_packages <- c("btools", "fantaxtic", "ampvis2", "tsnemicrobiota") # nombre del paquete

# Install CRAN packages
.inst <- cran_packages %in% installed.packages()
if(any(!.inst)) {
  install.packages(cran_packages[!.inst])
}
# Install Bioconductor packages
if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
.inst <- bioc_packages %in% installed.packages()
if(any(!.inst)) {
  BiocManager::install(bioc_packages[!.inst])
}
## Bioconductor version 3.12 (BiocManager 1.30.16), R 4.0.3 (2020-10-10)
## Installing package(s) 'philr'
## also installing the dependency 'ggtree'
## Warning in .inet_warning(msg): installation of package 'ggtree' had non-zero
## exit status
## Warning in .inet_warning(msg): installation of package 'philr' had non-zero exit
## status
## Updating HTML index of packages in '.Library'
## Making 'packages.html' ... done
## Old packages: 'aplot', 'class', 'fansi', 'fftwtools', 'foreign', 'MASS',
##   'nlme', 'nnet', 'Rcpp', 'reticulate', 'spatial', 'tidytree'
# Install GitHub packages
.inst <- git_source %in% installed.packages()
if(any(!.inst)) {
  devtools::install_github(git_source[!.inst])
}
## Skipping install of 'btools' from a github remote, the SHA1 (fa21d4ca) has not changed since last install.
##   Use `force = TRUE` to force installation
## Skipping install of 'fantaxtic' from a github remote, the SHA1 (3342beee) has not changed since last install.
##   Use `force = TRUE` to force installation
## Skipping install of 'ampvis2' from a github remote, the SHA1 (f48b708f) has not changed since last install.
##   Use `force = TRUE` to force installation
## Skipping install of 'tsnemicrobiota' from a github remote, the SHA1 (cffbab72) has not changed since last install.
##   Use `force = TRUE` to force installation
# Loading required packages
sapply(c(cran_packages, bioc_packages, git_packages), require, character.only = TRUE)
## Loading required package: bookdown
## Loading required package: knitr
## Loading required package: tidyverse
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──
## ✓ ggplot2 3.3.5     ✓ purrr   0.3.4
## ✓ tibble  3.1.6     ✓ dplyr   1.0.7
## ✓ tidyr   1.1.4     ✓ stringr 1.4.0
## ✓ readr   2.1.1     ✓ forcats 0.5.1
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
## Loading required package: plyr
## ------------------------------------------------------------------------------
## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)
## ------------------------------------------------------------------------------
## 
## Attaching package: 'plyr'
## 
## The following objects are masked from 'package:dplyr':
## 
##     arrange, count, desc, failwith, id, mutate, rename, summarise,
##     summarize
## 
## The following object is masked from 'package:purrr':
## 
##     compact
## 
## Loading required package: grid
## Loading required package: gridExtra
## 
## Attaching package: 'gridExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     combine
## 
## Loading required package: kableExtra
## 
## Attaching package: 'kableExtra'
## 
## The following object is masked from 'package:dplyr':
## 
##     group_rows
## 
## Loading required package: xtable
## Loading required package: ggpubr
## 
## Attaching package: 'ggpubr'
## 
## The following object is masked from 'package:plyr':
## 
##     mutate
## 
## Loading required package: dada2
## Loading required package: Rcpp
## Loading required package: DECIPHER
## Loading required package: Biostrings
## Loading required package: BiocGenerics
## Loading required package: parallel
## 
## Attaching package: 'BiocGenerics'
## 
## The following objects are masked from 'package:parallel':
## 
##     clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
##     clusterExport, clusterMap, parApply, parCapply, parLapply,
##     parLapplyLB, parRapply, parSapply, parSapplyLB
## 
## The following object is masked from 'package:gridExtra':
## 
##     combine
## 
## The following objects are masked from 'package:dplyr':
## 
##     combine, intersect, setdiff, union
## 
## The following objects are masked from 'package:stats':
## 
##     IQR, mad, sd, var, xtabs
## 
## The following objects are masked from 'package:base':
## 
##     anyDuplicated, append, as.data.frame, basename, cbind, colnames,
##     dirname, do.call, duplicated, eval, evalq, Filter, Find, get, grep,
##     grepl, intersect, is.unsorted, lapply, Map, mapply, match, mget,
##     order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
##     rbind, Reduce, rownames, sapply, setdiff, sort, table, tapply,
##     union, unique, unsplit, which.max, which.min
## 
## Loading required package: S4Vectors
## Loading required package: stats4
## 
## Attaching package: 'S4Vectors'
## 
## The following object is masked from 'package:plyr':
## 
##     rename
## 
## The following objects are masked from 'package:dplyr':
## 
##     first, rename
## 
## The following object is masked from 'package:tidyr':
## 
##     expand
## 
## The following object is masked from 'package:base':
## 
##     expand.grid
## 
## Loading required package: IRanges
## 
## Attaching package: 'IRanges'
## 
## The following object is masked from 'package:plyr':
## 
##     desc
## 
## The following objects are masked from 'package:dplyr':
## 
##     collapse, desc, slice
## 
## The following object is masked from 'package:purrr':
## 
##     reduce
## 
## The following object is masked from 'package:phyloseq':
## 
##     distance
## 
## Loading required package: XVector
## 
## Attaching package: 'XVector'
## 
## The following object is masked from 'package:plyr':
## 
##     compact
## 
## The following object is masked from 'package:purrr':
## 
##     compact
## 
## 
## Attaching package: 'Biostrings'
## 
## The following object is masked from 'package:base':
## 
##     strsplit
## 
## Loading required package: RSQLite
## Loading required package: phangorn
## Warning: package 'phangorn' was built under R version 4.0.5
## Loading required package: ape
## 
## Attaching package: 'ape'
## 
## The following object is masked from 'package:Biostrings':
## 
##     complement
## 
## The following object is masked from 'package:ggpubr':
## 
##     rotate
## 
## Loading required package: BiocManager
## Bioconductor version '3.12' is out-of-date; the current release version '3.14'
##   is available with R version '4.1'; see https://bioconductor.org/install
## Loading required package: DESeq2
## Loading required package: GenomicRanges
## Loading required package: GenomeInfoDb
## Loading required package: SummarizedExperiment
## Loading required package: MatrixGenerics
## Loading required package: matrixStats
## 
## Attaching package: 'matrixStats'
## 
## The following object is masked from 'package:plyr':
## 
##     count
## 
## The following object is masked from 'package:dplyr':
## 
##     count
## 
## 
## Attaching package: 'MatrixGenerics'
## 
## The following objects are masked from 'package:matrixStats':
## 
##     colAlls, colAnyNAs, colAnys, colAvgsPerRowSet, colCollapse,
##     colCounts, colCummaxs, colCummins, colCumprods, colCumsums,
##     colDiffs, colIQRDiffs, colIQRs, colLogSumExps, colMadDiffs,
##     colMads, colMaxs, colMeans2, colMedians, colMins, colOrderStats,
##     colProds, colQuantiles, colRanges, colRanks, colSdDiffs, colSds,
##     colSums2, colTabulates, colVarDiffs, colVars, colWeightedMads,
##     colWeightedMeans, colWeightedMedians, colWeightedSds,
##     colWeightedVars, rowAlls, rowAnyNAs, rowAnys, rowAvgsPerColSet,
##     rowCollapse, rowCounts, rowCummaxs, rowCummins, rowCumprods,
##     rowCumsums, rowDiffs, rowIQRDiffs, rowIQRs, rowLogSumExps,
##     rowMadDiffs, rowMads, rowMaxs, rowMeans2, rowMedians, rowMins,
##     rowOrderStats, rowProds, rowQuantiles, rowRanges, rowRanks,
##     rowSdDiffs, rowSds, rowSums2, rowTabulates, rowVarDiffs, rowVars,
##     rowWeightedMads, rowWeightedMeans, rowWeightedMedians,
##     rowWeightedSds, rowWeightedVars
## 
## Loading required package: Biobase
## Welcome to Bioconductor
## 
##     Vignettes contain introductory material; view with
##     'browseVignettes()'. To cite Bioconductor, see
##     'citation("Biobase")', and for packages 'citation("pkgname")'.
## 
## 
## Attaching package: 'Biobase'
## 
## The following object is masked from 'package:MatrixGenerics':
## 
##     rowMedians
## 
## The following objects are masked from 'package:matrixStats':
## 
##     anyMissing, rowMedians
## 
## The following object is masked from 'package:phyloseq':
## 
##     sampleNames
## 
## Loading required package: microbiome
## 
## microbiome R package (microbiome.github.com)
##     
## 
## 
##  Copyright (C) 2011-2020 Leo Lahti, 
##     Sudarshan Shetty et al. <microbiome.github.io>
## 
## 
## Attaching package: 'microbiome'
## 
## The following object is masked from 'package:SummarizedExperiment':
## 
##     coverage
## 
## The following object is masked from 'package:GenomicRanges':
## 
##     coverage
## 
## The following object is masked from 'package:phangorn':
## 
##     diversity
## 
## The following object is masked from 'package:Biostrings':
## 
##     coverage
## 
## The following objects are masked from 'package:IRanges':
## 
##     coverage, transform
## 
## The following object is masked from 'package:S4Vectors':
## 
##     transform
## 
## The following object is masked from 'package:ggplot2':
## 
##     alpha
## 
## The following object is masked from 'package:base':
## 
##     transform
## 
## Loading required package: philr
## Warning in library(package, lib.loc = lib.loc, character.only = TRUE,
## logical.return = TRUE, : there is no package called 'philr'
## Loading required package: btools
## Loading required package: fantaxtic
## Loading required package: reshape2
## 
## Attaching package: 'reshape2'
## 
## The following object is masked from 'package:tidyr':
## 
##     smiths
## 
## Loading required package: ampvis2
## Loading required package: tsnemicrobiota
##       bookdown          knitr      tidyverse           plyr           grid 
##           TRUE           TRUE           TRUE           TRUE           TRUE 
##      gridExtra     kableExtra         xtable         ggpubr       phyloseq 
##           TRUE           TRUE           TRUE           TRUE           TRUE 
##          dada2       DECIPHER       phangorn         ggpubr    BiocManager 
##           TRUE           TRUE           TRUE           TRUE           TRUE 
##         DESeq2     microbiome          philr         btools      fantaxtic 
##           TRUE           TRUE          FALSE           TRUE           TRUE 
##        ampvis2 tsnemicrobiota 
##           TRUE           TRUE

#Extract sample names

miseq_path <- "./" # CHANGE to the directory containing the fastq files after unzipping.
# Sort ensures forward/reverse reads are in same order
fnFs <- sort(list.files(miseq_path, pattern="_R1_001.fastq"))

# Extract sample names, assuming filenames have format: SAMPLENAME_XXX.fastq
sampleNames <- sapply(strsplit(fnFs, "_"), `[`, 1)
# Specify the full path to the fnFs
fnFs <- file.path(miseq_path, fnFs)

fnFs[1:3]
## [1] ".//CUN02-V2_1_L001_R1_001.fastq" ".//CUN02-V3_4_L001_R1_001.fastq"
## [3] ".//CUN03-V2_7_L001_R1_001.fastq"
#plot quality profile two samples
plotQualityProfile(fnFs[1:2])
## Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> =
## "none")` instead.

#plot quality profile 
plotQualityProfile(fnFs[c(1,11)])
## Warning: `guides(<scale> = FALSE)` is deprecated. Please use `guides(<scale> =
## "none")` instead.

#names(filtFs) <- sample.names
sample.names <- sapply(strsplit(basename(fnFs), "_"), `[`, 1)
head(sample.names)
## [1] "CUN02-V2" "CUN02-V3" "CUN03-V2" "CUN03-V3" "HCB01-V2" "HCB01-V3"
filt_path <- file.path(miseq_path, "filtered") # Place filtered files in filtered/ subdirectory
if(!file_test("-d", filt_path)) dir.create(filt_path)
filtFs <- file.path(filt_path, paste0(sampleNames, "_F_filt.fastq.gz"))

#Fiter and trim

outF1s<- filterAndTrim(fnFs, filtFs, truncLen=288, maxN=0, maxEE=2, truncQ=2, rm.phix=FALSE, compress=FALSE, verbose=TRUE)
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/CUN02-V2_F_filt.fastq.gz
## Read in 77582, output 55882 (72%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/CUN02-V3_F_filt.fastq.gz
## Read in 78500, output 56171 (71.6%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/CUN03-V2_F_filt.fastq.gz
## Read in 73044, output 50919 (69.7%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/CUN03-V3_F_filt.fastq.gz
## Read in 85148, output 60267 (70.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB01-V2_F_filt.fastq.gz
## Read in 77942, output 54420 (69.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB01-V3_F_filt.fastq.gz
## Read in 76092, output 51823 (68.1%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB02-V2_F_filt.fastq.gz
## Read in 62110, output 43706 (70.4%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB02-V3_F_filt.fastq.gz
## Read in 79907, output 56642 (70.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB03-V2_F_filt.fastq.gz
## Read in 90718, output 59814 (65.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB03-V3_F_filt.fastq.gz
## Read in 74075, output 49537 (66.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB04-V2_F_filt.fastq.gz
## Read in 78515, output 54934 (70%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB04-V3_F_filt.fastq.gz
## Read in 87173, output 60951 (69.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB06-V2_F_filt.fastq.gz
## Read in 86527, output 60214 (69.6%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB08-V2_F_filt.fastq.gz
## Read in 73143, output 49951 (68.3%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB09-V2_F_filt.fastq.gz
## Read in 92119, output 61560 (66.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB10-V2_F_filt.fastq.gz
## Read in 50022, output 35458 (70.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB12-V2_F_filt.fastq.gz
## Read in 73294, output 50478 (68.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HCB13-V3_F_filt.fastq.gz
## Read in 80869, output 55595 (68.7%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA01-V2_F_filt.fastq.gz
## Read in 67690, output 46666 (68.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA01-V3_F_filt.fastq.gz
## Read in 81275, output 57916 (71.3%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA02-V2_F_filt.fastq.gz
## Read in 81425, output 55645 (68.3%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA02-V3_F_filt.fastq.gz
## Read in 77772, output 53549 (68.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA03-V2_F_filt.fastq.gz
## Read in 78250, output 53339 (68.2%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA03-V3_F_filt.fastq.gz
## Read in 86178, output 60153 (69.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA04-V2_F_filt.fastq.gz
## Read in 103598, output 71280 (68.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA04-V3_F_filt.fastq.gz
## Read in 74184, output 50179 (67.6%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVA05-V3_F_filt.fastq.gz
## Read in 94023, output 63596 (67.6%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVR01-V2_F_filt.fastq.gz
## Read in 48149, output 32090 (66.6%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/HVR02-V2_F_filt.fastq.gz
## Read in 66957, output 44071 (65.8%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS01-V2_F_filt.fastq.gz
## Read in 82565, output 57700 (69.9%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS03-V2_F_filt.fastq.gz
## Read in 244723, output 168013 (68.7%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS06-V2_F_filt.fastq.gz
## Read in 81278, output 55876 (68.7%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS08-V2_F_filt.fastq.gz
## Read in 44945, output 28090 (62.5%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS09-V2_F_filt.fastq.gz
## Read in 60124, output 38568 (64.1%) filtered sequences.
## Overwriting file:/home/eduvb/qiime2/Rstudio/filtered/ISBS10-V2_F_filt.fastq.gz
## Read in 72508, output 48676 (67.1%) filtered sequences.
head(outF1s)
##                               reads.in reads.out
## CUN02-V2_1_L001_R1_001.fastq     77582     55882
## CUN02-V3_4_L001_R1_001.fastq     78500     56171
## CUN03-V2_7_L001_R1_001.fastq     73044     50919
## CUN03-V3_9_L001_R1_001.fastq     85148     60267
## HCB01-V2_13_L001_R1_001.fastq    77942     54420
## HCB01-V3_16_L001_R1_001.fastq    76092     51823

#Dereplication

derepFs <- derepFastq(filtFs, verbose=TRUE)
## Dereplicating sequence entries in Fastq file: .//filtered/CUN02-V2_F_filt.fastq.gz
## Encountered 19350 unique sequences from 55882 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/CUN02-V3_F_filt.fastq.gz
## Encountered 20061 unique sequences from 56171 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/CUN03-V2_F_filt.fastq.gz
## Encountered 20630 unique sequences from 50919 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/CUN03-V3_F_filt.fastq.gz
## Encountered 24207 unique sequences from 60267 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB01-V2_F_filt.fastq.gz
## Encountered 22104 unique sequences from 54420 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB01-V3_F_filt.fastq.gz
## Encountered 19675 unique sequences from 51823 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB02-V2_F_filt.fastq.gz
## Encountered 14022 unique sequences from 43706 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB02-V3_F_filt.fastq.gz
## Encountered 18830 unique sequences from 56642 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB03-V2_F_filt.fastq.gz
## Encountered 20981 unique sequences from 59814 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB03-V3_F_filt.fastq.gz
## Encountered 17912 unique sequences from 49537 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB04-V2_F_filt.fastq.gz
## Encountered 17659 unique sequences from 54934 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB04-V3_F_filt.fastq.gz
## Encountered 19467 unique sequences from 60951 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB06-V2_F_filt.fastq.gz
## Encountered 23786 unique sequences from 60214 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB08-V2_F_filt.fastq.gz
## Encountered 17432 unique sequences from 49951 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB09-V2_F_filt.fastq.gz
## Encountered 24622 unique sequences from 61560 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB10-V2_F_filt.fastq.gz
## Encountered 12603 unique sequences from 35458 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB12-V2_F_filt.fastq.gz
## Encountered 17518 unique sequences from 50478 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HCB13-V3_F_filt.fastq.gz
## Encountered 22507 unique sequences from 55595 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA01-V2_F_filt.fastq.gz
## Encountered 18268 unique sequences from 46666 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA01-V3_F_filt.fastq.gz
## Encountered 21055 unique sequences from 57916 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA02-V2_F_filt.fastq.gz
## Encountered 21588 unique sequences from 55645 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA02-V3_F_filt.fastq.gz
## Encountered 21872 unique sequences from 53549 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA03-V2_F_filt.fastq.gz
## Encountered 18099 unique sequences from 53339 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA03-V3_F_filt.fastq.gz
## Encountered 17148 unique sequences from 60153 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA04-V2_F_filt.fastq.gz
## Encountered 31215 unique sequences from 71280 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA04-V3_F_filt.fastq.gz
## Encountered 20455 unique sequences from 50179 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVA05-V3_F_filt.fastq.gz
## Encountered 25355 unique sequences from 63596 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVR01-V2_F_filt.fastq.gz
## Encountered 13935 unique sequences from 32090 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/HVR02-V2_F_filt.fastq.gz
## Encountered 17785 unique sequences from 44071 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS01-V2_F_filt.fastq.gz
## Encountered 21785 unique sequences from 57700 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS03-V2_F_filt.fastq.gz
## Encountered 56826 unique sequences from 168013 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS06-V2_F_filt.fastq.gz
## Encountered 25194 unique sequences from 55876 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS08-V2_F_filt.fastq.gz
## Encountered 12163 unique sequences from 28090 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS09-V2_F_filt.fastq.gz
## Encountered 13497 unique sequences from 38568 total sequences read.
## Dereplicating sequence entries in Fastq file: .//filtered/ISBS10-V2_F_filt.fastq.gz
## Encountered 17999 unique sequences from 48676 total sequences read.
# Name the derep-class objects by the sample names
names(derepFs) <- sampleNames
#Learn errors
errF <- learnErrors(filtFs, multithread=TRUE)
## 107478144 total bases in 373188 reads from 7 samples will be used for learning the error rates.
plotErrors(errF)

#DADA2

dadaFs <- dada(derepFs, err=errF, multithread=TRUE)
## Sample 1 - 55882 reads in 19350 unique sequences.
## Sample 2 - 56171 reads in 20061 unique sequences.
## Sample 3 - 50919 reads in 20630 unique sequences.
## Sample 4 - 60267 reads in 24207 unique sequences.
## Sample 5 - 54420 reads in 22104 unique sequences.
## Sample 6 - 51823 reads in 19675 unique sequences.
## Sample 7 - 43706 reads in 14022 unique sequences.
## Sample 8 - 56642 reads in 18830 unique sequences.
## Sample 9 - 59814 reads in 20981 unique sequences.
## Sample 10 - 49537 reads in 17912 unique sequences.
## Sample 11 - 54934 reads in 17659 unique sequences.
## Sample 12 - 60951 reads in 19467 unique sequences.
## Sample 13 - 60214 reads in 23786 unique sequences.
## Sample 14 - 49951 reads in 17432 unique sequences.
## Sample 15 - 61560 reads in 24622 unique sequences.
## Sample 16 - 35458 reads in 12603 unique sequences.
## Sample 17 - 50478 reads in 17518 unique sequences.
## Sample 18 - 55595 reads in 22507 unique sequences.
## Sample 19 - 46666 reads in 18268 unique sequences.
## Sample 20 - 57916 reads in 21055 unique sequences.
## Sample 21 - 55645 reads in 21588 unique sequences.
## Sample 22 - 53549 reads in 21872 unique sequences.
## Sample 23 - 53339 reads in 18099 unique sequences.
## Sample 24 - 60153 reads in 17148 unique sequences.
## Sample 25 - 71280 reads in 31215 unique sequences.
## Sample 26 - 50179 reads in 20455 unique sequences.
## Sample 27 - 63596 reads in 25355 unique sequences.
## Sample 28 - 32090 reads in 13935 unique sequences.
## Sample 29 - 44071 reads in 17785 unique sequences.
## Sample 30 - 57700 reads in 21785 unique sequences.
## Sample 31 - 168013 reads in 56826 unique sequences.
## Sample 32 - 55876 reads in 25194 unique sequences.
## Sample 33 - 28090 reads in 12163 unique sequences.
## Sample 34 - 38568 reads in 13497 unique sequences.
## Sample 35 - 48676 reads in 17999 unique sequences.
#Inspect dada2 object
dadaFs[[1]]
## dada-class: object describing DADA2 denoising results
## 383 sequence variants were inferred from 19350 input unique sequences.
## Key parameters: OMEGA_A = 1e-40, OMEGA_C = 1e-40, BAND_SIZE = 16

##Construct sequence table and remove chimeras

library(dada2)
library(Rcpp)
#BiocManager::install('taRifx')
library(taRifx)
## 
## Attaching package: 'taRifx'
## The following object is masked from 'package:SummarizedExperiment':
## 
##     shift
## The following object is masked from 'package:GenomicRanges':
## 
##     shift
## The following object is masked from 'package:IRanges':
## 
##     shift
## The following object is masked from 'package:S4Vectors':
## 
##     first
## The following objects are masked from 'package:dplyr':
## 
##     between, distinct, first, last
## The following object is masked from 'package:purrr':
## 
##     rep_along
mergers <-merge.list(dadaFs, derepFs)

seqtabAll <- makeSequenceTable(mergers[!grepl("Mock", names(mergers))])
table(nchar(getSequences(seqtabAll)))
## 
##   288 
## 10484
#remove chimeric sequences 
seqtabNoC <- removeBimeraDenovo(seqtabAll)

#Assign taxonomy <<>>

#SILVA Reference Genetic Database

#wget \
#  -O "silva_nr99_v138.1_train_set.fa.gz" \
#  "https://zenodo.org/record/4587955/files/silva_nr99_v138.1_train_set.fa.gz?download=1"

#wget \
#  -O "silva_species_assignment_v138.1.fa.gz" \
#  "https://zenodo.org/record/4587955/files/silva_species_assignment_v138.1.fa.gz?download=1"

tt <- assignTaxonomy(seqtabNoC, "silva_nr99_v138.1_train_set.fa.gz")
#tt.plus <- addSpecies(tt, "silva_species_assignment_v138.1.fa.gz", verbose=TRUE)
#load("tt.plus.Rda")

#unname(head(tt.plus))
unname(head(tt))
##      [,1]       [,2]             [,3]                  [,4]              
## [1,] "Bacteria" "Firmicutes"     "Bacilli"             "Lactobacillales" 
## [2,] "Bacteria" "Firmicutes"     "Bacilli"             "Lactobacillales" 
## [3,] "Bacteria" "Firmicutes"     "Clostridia"          "Lachnospirales"  
## [4,] "Bacteria" "Firmicutes"     "Bacilli"             "Lactobacillales" 
## [5,] "Bacteria" "Proteobacteria" "Gammaproteobacteria" "Enterobacterales"
## [6,] "Bacteria" "Proteobacteria" "Gammaproteobacteria" "Enterobacterales"
##      [,5]                 [,6]                  
## [1,] "Streptococcaceae"   "Streptococcus"       
## [2,] "Streptococcaceae"   "Streptococcus"       
## [3,] "Lachnospiraceae"    "Blautia"             
## [4,] "Streptococcaceae"   "Streptococcus"       
## [5,] "Enterobacteriaceae" "Escherichia-Shigella"
## [6,] "Enterobacteriaceae" "Escherichia-Shigella"
taxTab <- tt

#Construct phylogenetic tree

#using the DECIPHER R package
seqs <- getSequences(seqtabNoC)
names(seqs) <- seqs # This propagates to the tip labels of the tree
alignment <- AlignSeqs(DNAStringSet(seqs), anchor=NA,verbose=FALSE)


#phangorn R package is then used to construct a phylogenetic tree

phangAlign <- phyDat(as(alignment, "matrix"), type="DNA")
dm <- dist.ml(phangAlign)

treeNJ <- NJ(dm) # Note, tip order != sequence order

fit = pml(treeNJ, data=phangAlign)
## negative edges length changed to 0!
fitGTR <- update(fit, k=4, inv=0.2)
fitGTR <- optim.pml(fitGTR, model="GTR", optInv=TRUE, optGamma=TRUE,
        rearrangement = "stochastic", control = pml.control(trace = 0))
detach("package:phangorn", unload=TRUE)

##Combine data into a phyloseq object

#Load metadata file
samdf <- read.csv("Metadata_OTUs_mod_1.csv",header=TRUE)
rownames(samdf) <- samdf$SampleID

keep.cols <- c("SampleID",  "Edad_a_la_inclusion", "Sexo",  "Respuesta_100",    "aloTPH_previo_si_no",  "Citogenética_alto_riesgo", "Afect_Extramedular_PET_screening_si_no",   "Dias_producción",  "CRS_si_no", "Grado_máx_CRS",   "HOSPITAL", "Visit")

samdf <- samdf[rownames(seqtabAll), keep.cols]

# the sample metadata, the sequence taxonomies, and the phylogenetic tree – can now be combined into a single object
ps <- phyloseq(otu_table(seqtabNoC, taxa_are_rows=FALSE), 
               sample_data(samdf), 
               tax_table(taxTab),phy_tree(fitGTR$tree))

ps
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 6149 taxa and 35 samples ]
## sample_data() Sample Data:       [ 35 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 6149 taxa by 6 taxonomic ranks ]
## phy_tree()    Phylogenetic Tree: [ 6149 tips and 6147 internal nodes ]
#Loading required libraries

library(tidyverse)
#BiocManager::install("curatedMetagenomicData")
library(curatedMetagenomicData)
## Loading required package: AnnotationHub
## Loading required package: BiocFileCache
## Loading required package: dbplyr
## 
## Attaching package: 'dbplyr'
## The following objects are masked from 'package:dplyr':
## 
##     ident, sql
## 
## Attaching package: 'AnnotationHub'
## The following object is masked from 'package:Biobase':
## 
##     cache
## Loading required package: ExperimentHub
##                       _           _
##   ___ _   _ _ __ __ _| |_ ___  __| |
##  / __| | | | '__/ _` | __/ _ \/ _` |
## | (__| |_| | | | (_| | ||  __/ (_| |
##  \___|\__,_|_|  \__,_|\__\___|\__,_|
##  __  __      _                                         _
## |  \/  | ___| |_ __ _  __ _  ___ _ __   ___  _ __ ___ (_) ___
## | |\/| |/ _ \ __/ _` |/ _` |/ _ \ '_ \ / _ \| '_ ` _ \| |/ __|
## | |  | |  __/ || (_| | (_| |  __/ | | | (_) | | | | | | | (__
## |_|  |_|\___|\__\__,_|\__, |\___|_| |_|\___/|_| |_| |_|_|\___|
##  ____        _        |___/
## |  _ \  __ _| |_ __ _
## | | | |/ _` | __/ _` |
## | |_| | (_| | || (_| |
## |____/ \__,_|\__\__,_|
library(phyloseq)
library(DESeq2)
#BiocManager::install("apeglm")
library(apeglm)
#BiocManager::install("scran")
#BiocManager::install("zinbwave")
library(zinbwave)
## Loading required package: SingleCellExperiment
library(scran)
library(cowplot)
## 
## Attaching package: 'cowplot'
## The following object is masked from 'package:ggpubr':
## 
##     get_legend
#BiocManager::install("VennDiagram")
library(VennDiagram)
## Loading required package: futile.logger
## 
## Attaching package: 'VennDiagram'
## The following object is masked from 'package:ape':
## 
##     rotate
## The following object is masked from 'package:ggpubr':
## 
##     rotate
#BiocManager::install("pscl")
library(pscl)
## Classes and Methods for R developed in the
## Political Science Computational Laboratory
## Department of Political Science
## Stanford University
## Simon Jackman
## hurdle and zeroinfl functions by Achim Zeileis
library(MASS)
## 
## Attaching package: 'MASS'
## The following object is masked from 'package:dplyr':
## 
##     select

#####Phyloseq object#####

#Taxonomic Filtering
# Show available ranks in the dataset
rank_names(ps)
## [1] "Kingdom" "Phylum"  "Class"   "Order"   "Family"  "Genus"
## Create table, number of features for each phyla
table(tax_table(ps)[, "Phylum"], exclude = NULL)
## 
##  Actinobacteriota      Bacteroidota  Campylobacterota     Cyanobacteria 
##               622               749                 6                 9 
##  Desulfobacterota     Euryarchaeota        Firmicutes    Fusobacteriota 
##                21                 8              4290                12 
##   Patescibacteria    Proteobacteria      Synergistota Verrucomicrobiota 
##                18               246                12                30 
##              <NA> 
##               126
#features were annotated with a Phylum of NA. These features are probably artifacts in a dataset like this, and should be removed.
ps <- subset_taxa(ps, !is.na(Phylum) & !Phylum %in% c("", "uncharacterized"))

# Compute prevalence of each feature, store as data.frame
prevdf = apply(X = otu_table(ps),
               MARGIN = ifelse(taxa_are_rows(ps), yes = 1, no = 2),
               FUN = function(x){sum(x > 0)})
# Add taxonomy and total read counts to this data.frame
prevdf = data.frame(Prevalence = prevdf,
                    TotalAbundance = taxa_sums(ps),
                    tax_table(ps))
#Create a new phyloseq-object
psv2 <- ps
#Define levels as factor
sample_data(psv2)$Respuesta_100 <- as.factor(sample_data(psv2)$Respuesta_100)
levels(sample_data(psv2)$Respuesta_100)
## [1] "PR"   "sCR"  "VGPR"
#rename <<VGPR>> en <<sCR>>
levels(sample_data(psv2)$Respuesta_100)[3] <- "sCR"
levels(sample_data(psv2)$Respuesta_100)
## [1] "PR"  "sCR"
sample_data(psv2)$Respuesta_100
##  [1] sCR sCR sCR sCR PR  PR  PR  PR  sCR sCR sCR sCR sCR PR  sCR sCR sCR sCR PR 
## [20] PR  sCR sCR sCR sCR sCR sCR PR  sCR sCR sCR sCR PR  sCR PR  sCR
## Levels: PR sCR
#Create a new phyloseq-object filtered by <<Visit==V2>>

psv2 <- subset_samples(psv2, Visit=="V2" )
sample_data(psv2)
##            SampleID Edad_a_la_inclusion   Sexo Respuesta_100
## CUN02-V2   CUN02-V2                  43  MUJER           sCR
## CUN03-V2   CUN03-V2                  45 HOMBRE           sCR
## HCB01-V2   HCB01-V2                  66 HOMBRE            PR
## HCB02-V2   HCB02-V2                  66 HOMBRE            PR
## HCB03-V2   HCB03-V2                  60 HOMBRE           sCR
## HCB04-V2   HCB04-V2                  55  MUJER           sCR
## HCB06-V2   HCB06-V2                  73  MUJER           sCR
## HCB08-V2   HCB08-V2                  62  MUJER            PR
## HCB09-V2   HCB09-V2                  53 HOMBRE           sCR
## HCB10-V2   HCB10-V2                  52 HOMBRE           sCR
## HCB12-V2   HCB12-V2                  64 HOMBRE           sCR
## HVA01-V2   HVA01-V2                  58 HOMBRE            PR
## HVA02-V2   HVA02-V2                  67 HOMBRE           sCR
## HVA03-V2   HVA03-V2                  61 HOMBRE           sCR
## HVA04-V2   HVA04-V2                  74  MUJER           sCR
## HVR01-V2   HVR01-V2                  65 HOMBRE           sCR
## HVR02-V2   HVR02-V2                  46 HOMBRE           sCR
## ISBS01-V2 ISBS01-V2                  56 HOMBRE           sCR
## ISBS03-V2 ISBS03-V2                  70  MUJER           sCR
## ISBS06-V2 ISBS06-V2                  60 HOMBRE            PR
## ISBS08-V2 ISBS08-V2                  53 HOMBRE           sCR
## ISBS09-V2 ISBS09-V2                  65  MUJER            PR
## ISBS10-V2 ISBS10-V2                  68 HOMBRE           sCR
##           aloTPH_previo_si_no Citogenética_alto_riesgo
## CUN02-V2                   NO                       NO
## CUN03-V2                   NO                       SI
## HCB01-V2                   NO                       NO
## HCB02-V2                   NO                       SI
## HCB03-V2                   NO                       SI
## HCB04-V2                   NO                       NO
## HCB06-V2                   NO                       NO
## HCB08-V2                   NO                       NO
## HCB09-V2                   NO                       NO
## HCB10-V2                   SI                       NO
## HCB12-V2                   NO                       SI
## HVA01-V2                   NO                       NO
## HVA02-V2                   NO                       NO
## HVA03-V2                   SI                       NO
## HVA04-V2                   NO                       NO
## HVR01-V2                                            NO
## HVR02-V2                                            NO
## ISBS01-V2                  NO                       SI
## ISBS03-V2                  NO                         
## ISBS06-V2                  NO                       NO
## ISBS08-V2                  NO                         
## ISBS09-V2                  NO                       NO
## ISBS10-V2                  NO                       SI
##           Afect_Extramedular_PET_screening_si_no Dias_producción CRS_si_no
## CUN02-V2                                      SI             120        SI
## CUN03-V2                                      NO             100        SI
## HCB01-V2                                      SI             110        SI
## HCB02-V2                                      NO             110        SI
## HCB03-V2                                      SI             110        SI
## HCB04-V2                                      NO             110        SI
## HCB06-V2                                      NO              90        SI
## HCB08-V2                                      SI             100        SI
## HCB09-V2                                      NO              90        SI
## HCB10-V2                                      SI             100        SI
## HCB12-V2                                      NO             110        NO
## HVA01-V2                                      NO             120        SI
## HVA02-V2                                      SI             100        SI
## HVA03-V2                                      SI             140        SI
## HVA04-V2                                      NO             100        SI
## HVR01-V2                                      SI             110        SI
## HVR02-V2                                      SI             120        SI
## ISBS01-V2                                     SI             120        SI
## ISBS03-V2                                     SI             100        SI
## ISBS06-V2                                     NO             130        SI
## ISBS08-V2                                     SI             110        SI
## ISBS09-V2                                     SI             100        SI
## ISBS10-V2                                     SI             140        NO
##           Grado_máx_CRS HOSPITAL Visit
## CUN02-V2             10      CUN    V2
## CUN03-V2             20      CUN    V2
## HCB01-V2             20      HCB    V2
## HCB02-V2             20      HCB    V2
## HCB03-V2             10      HCB    V2
## HCB04-V2             10      HCB    V2
## HCB06-V2             10      HCB    V2
## HCB08-V2             10      HCB    V2
## HCB09-V2             10      HCB    V2
## HCB10-V2             10      HCB    V2
## HCB12-V2             NA      HCB    V2
## HVA01-V2             20      HVA    V2
## HVA02-V2             20      HVA    V2
## HVA03-V2             10      HVA    V2
## HVA04-V2             10      HVA    V2
## HVR01-V2             10      HVR    V2
## HVR02-V2             10      HVR    V2
## ISBS01-V2            10     ISBS    V2
## ISBS03-V2            10     ISBS    V2
## ISBS06-V2            10     ISBS    V2
## ISBS08-V2            20     ISBS    V2
## ISBS09-V2            10     ISBS    V2
## ISBS10-V2            NA     ISBS    V2
#Review the levels of the variable  <<Respuesta_100>>
levels(sample_data(psv2)$Respuesta_100)
## [1] "PR"  "sCR"
#Summary variable <<Respuesta_100>>
summary(sample_data(psv2)$Respuesta_100)
##  PR sCR 
##   6  17
list(sample_data(psv2)$Visit)
## [[1]]
##  [1] "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2"
## [16] "V2" "V2" "V2" "V2" "V2" "V2" "V2" "V2"
#Filtered phyloseq object by patients visit <<V2>>
psv2
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 6023 taxa and 23 samples ]
## sample_data() Sample Data:       [ 23 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 6023 taxa by 6 taxonomic ranks ]
## phy_tree()    Phylogenetic Tree: [ 6023 tips and 6021 internal nodes ]

#Taxonomic prevalence

# Show available ranks in the dataset
rank_names(psv2)
## [1] "Kingdom" "Phylum"  "Class"   "Order"   "Family"  "Genus"
## Create table, number of features for each phyla
table(tax_table(psv2)[, "Phylum"], exclude = NULL)
## 
##  Actinobacteriota      Bacteroidota  Campylobacterota     Cyanobacteria 
##               622               749                 6                 9 
##  Desulfobacterota     Euryarchaeota        Firmicutes    Fusobacteriota 
##                21                 8              4290                12 
##   Patescibacteria    Proteobacteria      Synergistota Verrucomicrobiota 
##                18               246                12                30
#features were annotated with a Phylum of NA. These features are probably artifacts in a dataset like this, and should be removed.
psv2 <- subset_taxa(psv2, !is.na(Phylum) & !Phylum %in% c("", "uncharacterized"))


# compute prevalence for each feature and store it in a data frame
prevdf = apply(X = otu_table(psv2),
               MARGIN = ifelse(taxa_are_rows(psv2), yes = 1, no = 2),
               FUN = function(x){sum(x > 0)})
# add the taxonomy
prevdf = data.frame(Prevalence = prevdf,
                    TotalAbundance = taxa_sums(psv2),
                    tax_table(psv2))

plyr::ddply(prevdf, "Phylum", function(df1){cbind(mean(df1$Prevalence),sum(df1$Prevalence))}) -> dfprev
kable(dfprev)
Phylum 1 2
Actinobacteriota 1.1752412 731
Bacteroidota 1.2616822 945
Campylobacterota 0.8333333 5
Cyanobacteria 0.8888889 8
Desulfobacterota 1.1904762 25
Euryarchaeota 2.5000000 20
Firmicutes 1.2559441 5388
Fusobacteriota 0.5000000 6
Patescibacteria 0.8333333 15
Proteobacteria 1.1260163 277
Synergistota 1.0000000 12
Verrucomicrobiota 2.2666667 68
# Show available ranks in the dataset
rank_names(psv2)
## [1] "Kingdom" "Phylum"  "Class"   "Order"   "Family"  "Genus"
## Create table, number of features for each phyla
table(tax_table(psv2)[, "Phylum"], exclude = NULL)
## 
##  Actinobacteriota      Bacteroidota  Campylobacterota     Cyanobacteria 
##               622               749                 6                 9 
##  Desulfobacterota     Euryarchaeota        Firmicutes    Fusobacteriota 
##                21                 8              4290                12 
##   Patescibacteria    Proteobacteria      Synergistota Verrucomicrobiota 
##                18               246                12                30
#features were annotated with a Phylum of NA. These features are probably artifacts in a dataset like this, and should be removed.
psv2 <- subset_taxa(psv2, !is.na(Phylum) & !Phylum %in% c("", "uncharacterized"))


# Computamos prevalencia para cada feature y la guardamos en un data frame
prevdf = apply(X = otu_table(psv2),
               MARGIN = ifelse(taxa_are_rows(psv2), yes = 1, no = 2),
               FUN = function(x){sum(x > 0)})
# Le agregamos la taxonomía
prevdf = data.frame(Prevalence = prevdf,
                    TotalAbundance = taxa_sums(psv2),
                    tax_table(psv2))

plyr::ddply(prevdf, "Phylum", function(df1){cbind(mean(df1$Prevalence),sum(df1$Prevalence))}) -> dfprev
kable(dfprev)
Phylum 1 2
Actinobacteriota 1.1752412 731
Bacteroidota 1.2616822 945
Campylobacterota 0.8333333 5
Cyanobacteria 0.8888889 8
Desulfobacterota 1.1904762 25
Euryarchaeota 2.5000000 20
Firmicutes 1.2559441 5388
Fusobacteriota 0.5000000 6
Patescibacteria 0.8333333 15
Proteobacteria 1.1260163 277
Synergistota 1.0000000 12
Verrucomicrobiota 2.2666667 68
#La columna 1 representa la media de read counts para ese Phylum, mientras que la columna 2 representa la suma
# Define phyla to filter
filterPhyla = c("Campylobacterota", "Cyanobacteria", "Halanaerobiaeota", "Fusobacteriota", "Synergistota" )
psv2 = subset_taxa(psv2, !Phylum %in% filterPhyla)
psv2
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 5984 taxa and 23 samples ]
## sample_data() Sample Data:       [ 23 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 5984 taxa by 6 taxonomic ranks ]
## phy_tree()    Phylogenetic Tree: [ 5984 tips and 5982 internal nodes ]

#Filter taxa

# ilter taxa according to a threshold of mean number of _read counts_, in this case 1e-5
psd2 <- filter_taxa(psv2, function(x) mean(x) > 1e-5, TRUE)

# remove taxa that are not observed more than X times in at least 10% of the samples
psd3 <- filter_taxa(psd2, function(x) sum(x > 2) > (0.1*length(x)), TRUE)

# filter samples with less than 1000 reads
psd4 = prune_samples(sample_sums(psd3) > 1000, psd3)

psd4
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 557 taxa and 23 samples ]
## sample_data() Sample Data:       [ 23 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 557 taxa by 6 taxonomic ranks ]
## phy_tree()    Phylogenetic Tree: [ 557 tips and 555 internal nodes ]
# We select the taxa of interest (ALLS)
prevdf1 = subset(prevdf, Phylum %in% get_taxa_unique(psd4, "Phylum"))
ggplot(prevdf1, aes(TotalAbundance, Prevalence / nsamples(psv2),color=Phylum)) +
  # Agregamos una línea para nuestro umbral
  geom_hline(yintercept = 0.05, alpha = 0.5, linetype = 2) +  geom_point(size = 2, alpha = 0.7) +
  scale_x_log10() +  xlab("Total Abundance") + ylab("Prevalence [Frac. Samples]") +
  facet_wrap(~Phylum) + theme(legend.position="none")
## Warning: Transformation introduced infinite values in continuous x-axis

# replace the sequences with a generic name
psd5 <- psd4
taxa_names(psd5) <- paste0("ASV", seq(ntaxa(psd5)))
head(taxa_names(psd5))
## [1] "ASV1" "ASV2" "ASV3" "ASV4" "ASV5" "ASV6"
#plot the distribution of read counts by sample

sample_sum_df <- data.frame(sum = sample_sums(psd5))

ggplot(sample_sum_df, aes(x = sum)) + 
  geom_histogram(color = "black", fill = "grey", binwidth = 2500) +
  ggtitle("Distribution of sample sequencing depth") + 
  xlab("Read counts") +
  theme(axis.title.y = element_blank()) 

#calculate rarefaction curves for each sample

# First  load  scripts 
scripts <- c("graphical_methods.R",
             "tree_methods.R",
             "plot_merged_trees.R",
             "specificity_methods.R",
             "ternary_plot.R",
             "richness.R",
             "edgePCA.R",
             "copy_number_correction.R",
             "import_frogs.R",
             "prevalence.R",
             "compute_niche.R")
urls <- paste0("https://raw.githubusercontent.com/mahendra-mariadassou/phyloseq-extended/master/R/", scripts)

for (url in urls) {
  source(url)
}
## Loading required package: scales
## 
## Attaching package: 'scales'
## The following object is masked from 'package:microbiome':
## 
##     alpha
## The following object is masked from 'package:purrr':
## 
##     discard
## The following object is masked from 'package:readr':
## 
##     col_factor
## Loading required package: vegan
## Loading required package: permute
## Warning: package 'permute' was built under R version 4.0.5
## Loading required package: lattice
## This is vegan 2.5-7
## 
## Attaching package: 'vegan'
## The following object is masked from 'package:microbiome':
## 
##     diversity
#Plot the rarefaction curves

p <- ggrare(psd5, step = 100, color = "Respuesta_100", label = "SampleID", se = TRUE)
## rarefying sample CUN02-V2
## rarefying sample CUN03-V2
## rarefying sample HCB01-V2
## rarefying sample HCB02-V2
## rarefying sample HCB03-V2
## rarefying sample HCB04-V2
## rarefying sample HCB06-V2
## rarefying sample HCB08-V2
## rarefying sample HCB09-V2
## rarefying sample HCB10-V2
## rarefying sample HCB12-V2
## rarefying sample HVA01-V2
## rarefying sample HVA02-V2
## rarefying sample HVA03-V2
## rarefying sample HVA04-V2
## rarefying sample HVR01-V2
## rarefying sample HVR02-V2
## rarefying sample ISBS01-V2
## rarefying sample ISBS03-V2
## rarefying sample ISBS06-V2
## rarefying sample ISBS08-V2
## rarefying sample ISBS09-V2
## rarefying sample ISBS10-V2

(p <- p + facet_wrap(~Respuesta_100))

Table of read counts

#For visualization in HTML format
library("kableExtra")

head(otu_table(psd5)) %>%
  kable(format = "html", col.names = colnames(otu_table(psd5))) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "350px")
ASV1 ASV2 ASV3 ASV4 ASV5 ASV6 ASV7 ASV8 ASV9 ASV10 ASV11 ASV12 ASV13 ASV14 ASV15 ASV16 ASV17 ASV18 ASV19 ASV20 ASV21 ASV22 ASV23 ASV24 ASV25 ASV26 ASV27 ASV28 ASV29 ASV30 ASV31 ASV32 ASV33 ASV34 ASV35 ASV36 ASV37 ASV38 ASV39 ASV40 ASV41 ASV42 ASV43 ASV44 ASV45 ASV46 ASV47 ASV48 ASV49 ASV50 ASV51 ASV52 ASV53 ASV54 ASV55 ASV56 ASV57 ASV58 ASV59 ASV60 ASV61 ASV62 ASV63 ASV64 ASV65 ASV66 ASV67 ASV68 ASV69 ASV70 ASV71 ASV72 ASV73 ASV74 ASV75 ASV76 ASV77 ASV78 ASV79 ASV80 ASV81 ASV82 ASV83 ASV84 ASV85 ASV86 ASV87 ASV88 ASV89 ASV90 ASV91 ASV92 ASV93 ASV94 ASV95 ASV96 ASV97 ASV98 ASV99 ASV100 ASV101 ASV102 ASV103 ASV104 ASV105 ASV106 ASV107 ASV108 ASV109 ASV110 ASV111 ASV112 ASV113 ASV114 ASV115 ASV116 ASV117 ASV118 ASV119 ASV120 ASV121 ASV122 ASV123 ASV124 ASV125 ASV126 ASV127 ASV128 ASV129 ASV130 ASV131 ASV132 ASV133 ASV134 ASV135 ASV136 ASV137 ASV138 ASV139 ASV140 ASV141 ASV142 ASV143 ASV144 ASV145 ASV146 ASV147 ASV148 ASV149 ASV150 ASV151 ASV152 ASV153 ASV154 ASV155 ASV156 ASV157 ASV158 ASV159 ASV160 ASV161 ASV162 ASV163 ASV164 ASV165 ASV166 ASV167 ASV168 ASV169 ASV170 ASV171 ASV172 ASV173 ASV174 ASV175 ASV176 ASV177 ASV178 ASV179 ASV180 ASV181 ASV182 ASV183 ASV184 ASV185 ASV186 ASV187 ASV188 ASV189 ASV190 ASV191 ASV192 ASV193 ASV194 ASV195 ASV196 ASV197 ASV198 ASV199 ASV200 ASV201 ASV202 ASV203 ASV204 ASV205 ASV206 ASV207 ASV208 ASV209 ASV210 ASV211 ASV212 ASV213 ASV214 ASV215 ASV216 ASV217 ASV218 ASV219 ASV220 ASV221 ASV222 ASV223 ASV224 ASV225 ASV226 ASV227 ASV228 ASV229 ASV230 ASV231 ASV232 ASV233 ASV234 ASV235 ASV236 ASV237 ASV238 ASV239 ASV240 ASV241 ASV242 ASV243 ASV244 ASV245 ASV246 ASV247 ASV248 ASV249 ASV250 ASV251 ASV252 ASV253 ASV254 ASV255 ASV256 ASV257 ASV258 ASV259 ASV260 ASV261 ASV262 ASV263 ASV264 ASV265 ASV266 ASV267 ASV268 ASV269 ASV270 ASV271 ASV272 ASV273 ASV274 ASV275 ASV276 ASV277 ASV278 ASV279 ASV280 ASV281 ASV282 ASV283 ASV284 ASV285 ASV286 ASV287 ASV288 ASV289 ASV290 ASV291 ASV292 ASV293 ASV294 ASV295 ASV296 ASV297 ASV298 ASV299 ASV300 ASV301 ASV302 ASV303 ASV304 ASV305 ASV306 ASV307 ASV308 ASV309 ASV310 ASV311 ASV312 ASV313 ASV314 ASV315 ASV316 ASV317 ASV318 ASV319 ASV320 ASV321 ASV322 ASV323 ASV324 ASV325 ASV326 ASV327 ASV328 ASV329 ASV330 ASV331 ASV332 ASV333 ASV334 ASV335 ASV336 ASV337 ASV338 ASV339 ASV340 ASV341 ASV342 ASV343 ASV344 ASV345 ASV346 ASV347 ASV348 ASV349 ASV350 ASV351 ASV352 ASV353 ASV354 ASV355 ASV356 ASV357 ASV358 ASV359 ASV360 ASV361 ASV362 ASV363 ASV364 ASV365 ASV366 ASV367 ASV368 ASV369 ASV370 ASV371 ASV372 ASV373 ASV374 ASV375 ASV376 ASV377 ASV378 ASV379 ASV380 ASV381 ASV382 ASV383 ASV384 ASV385 ASV386 ASV387 ASV388 ASV389 ASV390 ASV391 ASV392 ASV393 ASV394 ASV395 ASV396 ASV397 ASV398 ASV399 ASV400 ASV401 ASV402 ASV403 ASV404 ASV405 ASV406 ASV407 ASV408 ASV409 ASV410 ASV411 ASV412 ASV413 ASV414 ASV415 ASV416 ASV417 ASV418 ASV419 ASV420 ASV421 ASV422 ASV423 ASV424 ASV425 ASV426 ASV427 ASV428 ASV429 ASV430 ASV431 ASV432 ASV433 ASV434 ASV435 ASV436 ASV437 ASV438 ASV439 ASV440 ASV441 ASV442 ASV443 ASV444 ASV445 ASV446 ASV447 ASV448 ASV449 ASV450 ASV451 ASV452 ASV453 ASV454 ASV455 ASV456 ASV457 ASV458 ASV459 ASV460 ASV461 ASV462 ASV463 ASV464 ASV465 ASV466 ASV467 ASV468 ASV469 ASV470 ASV471 ASV472 ASV473 ASV474 ASV475 ASV476 ASV477 ASV478 ASV479 ASV480 ASV481 ASV482 ASV483 ASV484 ASV485 ASV486 ASV487 ASV488 ASV489 ASV490 ASV491 ASV492 ASV493 ASV494 ASV495 ASV496 ASV497 ASV498 ASV499 ASV500 ASV501 ASV502 ASV503 ASV504 ASV505 ASV506 ASV507 ASV508 ASV509 ASV510 ASV511 ASV512 ASV513 ASV514 ASV515 ASV516 ASV517 ASV518 ASV519 ASV520 ASV521 ASV522 ASV523 ASV524 ASV525 ASV526 ASV527 ASV528 ASV529 ASV530 ASV531 ASV532 ASV533 ASV534 ASV535 ASV536 ASV537 ASV538 ASV539 ASV540 ASV541 ASV542 ASV543 ASV544 ASV545 ASV546 ASV547 ASV548 ASV549 ASV550 ASV551 ASV552 ASV553 ASV554 ASV555 ASV556 ASV557
CUN02-V2 0 0 528 0 0 0 447 1933 0 0 1863 272 0 213 0 0 1121 0 0 0 747 0 0 0 0 178 0 0 0 0 0 0 0 0 0 0 0 617 0 0 0 0 771 297 547 0 0 220 0 158 0 0 0 0 0 0 1099 718 334 110 188 0 0 0 0 0 745 0 0 0 0 0 0 0 0 0 0 0 380 0 0 0 0 0 0 0 1020 397 0 120 0 264 0 128 135 0 0 0 0 0 1983 2081 0 127 0 449 0 0 0 0 0 490 0 0 0 0 0 0 0 121 445 106 643 0 0 0 0 85 0 484 0 0 0 375 12 0 0 0 0 0 0 0 0 303 0 0 117 79 381 0 0 0 0 0 42 75 0 143 0 459 1247 0 0 357 388 68 0 0 0 0 225 0 395 0 0 0 0 0 0 245 0 0 135 0 0 0 341 79 0 0 0 57 0 0 0 0 0 221 0 0 185 852 215 0 93 0 157 0 0 0 0 0 0 0 59 0 0 0 0 0 0 52 0 0 0 0 0 204 0 90 0 0 44 0 118 217 0 0 0 0 0 104 0 0 0 28 0 0 0 0 0 0 0 0 0 88 0 0 91 0 0 0 0 0 0 0 283 0 0 0 0 0 0 0 0 0 56 0 0 0 0 0 125 108 0 0 0 104 0 269 0 0 0 0 0 0 0 0 0 0 0 0 84 0 0 267 0 69 0 49 591 111 0 41 0 0 0 0 0 0 89 0 0 0 0 134 0 53 0 0 0 0 0 0 51 0 141 0 0 0 0 0 0 0 0 0 0 0 0 0 274 104 125 0 0 0 0 0 0 181 0 0 0 0 0 0 0 44 0 176 171 125 0 0 47 0 0 112 0 0 0 0 0 0 0 0 0 0 103 0 0 0 0 0 0 0 0 0 0 0 42 0 0 60 0 0 0 0 0 0 0 0 0 0 0 0 0 144 0 0 0 0 0 149 0 0 0 0 0 0 0 87 0 0 108 0 0 0 0 0 0 0 0 52 0 0 0 0 0 0 0 0 0 0 6 0 0 0 0 61 168 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 9 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 59 0 0 0 0 0 46 0 0 87 0 0 0 0 59 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
CUN03-V2 243 241 183 183 481 454 159 0 141 308 0 142 737 86 230 170 0 110 710 133 0 100 492 104 125 170 189 1340 0 312 432 83 76 1391 111 69 327 78 23 0 0 90 0 141 44 0 0 0 65 0 205 130 924 0 0 205 94 141 0 0 0 0 0 0 95 212 135 13 0 0 0 158 0 184 152 0 0 0 0 764 60 186 0 197 0 184 97 0 0 0 0 29 0 110 0 117 64 0 0 0 0 0 166 51 0 91 0 0 0 0 0 110 139 100 0 0 69 0 60 0 0 0 82 50 114 569 0 0 182 130 0 150 91 86 0 124 33 0 0 0 60 107 0 91 57 0 137 0 0 0 73 0 168 0 33 0 88 42 0 50 0 0 56 186 227 0 0 0 0 0 177 0 269 0 0 0 146 156 0 0 0 0 87 130 0 0 207 0 51 0 0 0 0 106 76 0 370 216 0 0 0 0 113 302 91 0 0 0 12 0 0 0 0 0 108 0 60 0 0 0 0 0 0 161 0 0 154 197 0 0 0 0 0 0 0 134 61 0 0 0 0 0 0 113 0 90 0 0 123 0 0 0 0 0 0 96 139 0 0 0 157 0 0 0 0 0 0 0 202 0 0 0 0 0 0 0 59 0 0 0 0 0 104 0 58 0 46 108 0 0 83 0 160 0 123 0 0 0 0 0 0 0 0 102 0 0 0 0 38 0 0 0 0 0 0 0 0 97 13 0 0 113 68 86 0 0 0 74 0 0 0 0 145 0 106 0 0 0 0 115 0 57 0 0 0 0 83 58 0 83 0 0 0 0 85 0 12 0 84 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 61 0 0 0 180 0 0 59 0 0 0 57 0 0 0 0 178 0 146 0 0 0 0 0 0 0 0 187 141 0 0 83 0 27 0 0 0 0 0 0 0 0 0 0 0 0 52 0 0 0 32 0 0 37 149 0 0 0 97 0 0 93 0 114 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 119 0 51 0 0 0 0 0 0 0 0 0 0 0 0 61 0 0 94 0 0 0 0 0 0 33 123 0 0 52 0 0 0 39 0 0 0 0 0 0 0 0 0 0 71 0 0 0 0 0 83 152 0 0 0 158 0 0 0 7 0 0 0 0 0 0 17 0 0 54 0 0 3 0 51 60 0 0 39 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 34 0 7 0 0 0
HCB01-V2 131 136 408 119 0 0 372 0 59 0 0 271 0 220 0 300 0 58 0 352 0 363 0 215 310 247 0 1192 0 1145 0 170 0 1185 198 0 1078 0 0 0 0 535 0 223 0 0 0 0 193 182 774 498 854 0 0 229 279 107 0 88 0 0 638 249 0 248 90 617 187 0 55 83 79 221 205 0 0 0 0 706 406 154 0 634 0 218 225 0 0 214 0 0 0 117 0 139 296 0 223 0 0 0 114 220 0 66 0 0 0 0 0 0 87 0 0 376 141 0 0 144 89 168 179 0 0 449 0 0 114 0 190 149 140 63 122 444 110 0 0 0 0 81 0 0 0 0 0 139 0 0 106 0 0 0 62 71 82 175 155 139 0 0 0 0 0 115 89 0 141 130 0 0 0 0 0 0 107 78 0 0 0 0 0 0 0 90 0 0 0 79 0 114 0 0 0 101 303 0 0 0 82 0 0 310 0 316 172 0 0 0 64 0 0 304 0 0 41 0 0 0 294 221 0 0 101 0 0 0 0 223 0 0 0 0 0 0 0 76 77 0 221 0 284 0 0 0 0 0 0 337 224 140 0 38 313 65 0 96 0 0 0 0 204 60 0 29 0 0 217 170 0 0 142 0 115 0 0 97 0 0 153 68 0 48 0 0 218 0 97 0 0 159 0 117 0 166 0 0 0 107 0 0 0 273 137 0 63 0 0 0 0 0 142 136 218 196 98 0 160 0 0 0 0 108 86 0 0 0 42 0 33 141 0 0 0 0 0 151 74 0 54 0 0 0 0 0 0 0 60 0 0 0 0 0 0 267 0 0 0 0 73 95 27 0 0 0 36 0 0 0 0 0 0 0 0 0 102 0 0 0 0 0 0 117 78 0 0 0 127 0 0 0 0 87 87 227 0 0 0 43 0 0 0 0 0 0 109 0 97 0 0 0 0 0 0 57 0 0 0 0 0 0 0 0 0 100 0 0 0 0 0 0 0 0 0 0 0 62 0 30 0 0 0 0 0 0 0 0 0 0 0 0 50 0 0 0 85 0 0 0 0 0 0 0 0 0 35 0 0 0 71 0 0 0 0 0 0 0 95 0 0 0 0 0 0 0 0 0 0 0 0 102 31 0 110 0 0 0 0 43 0 0 0 0 0 0 0 41 0 0 0 0 0 0 0 0 0 61 35 35 0 71 0 0 0 0 0 0 67 0 0 75 0 81 37 0 25 10 0 0 0 0 0 0 0 15 13 0 0 0 0 0 0 0 6 0 6
HCB02-V2 250 0 0 0 0 0 0 413 54 0 396 62 0 0 0 0 277 0 0 0 244 0 0 0 0 0 0 0 0 124 0 0 0 0 39 0 0 0 0 0 0 0 177 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 273 308 548 114 0 0 0 0 0 112 0 0 563 0 28 0 0 0 106 86 191 0 0 0 0 0 0 0 0 0 0 0 0 0 196 0 0 0 0 207 0 0 0 349 0 0 0 0 0 0 0 62 100 0 0 0 0 0 0 0 0 0 0 0 0 135 0 0 127 0 0 0 259 0 0 0 0 0 101 64 0 0 0 0 0 0 239 0 0 0 0 0 0 0 457 0 0 0 0 80 0 0 249 0 0 0 0 0 0 1639 0 0 0 0 0 0 0 427 0 0 0 0 0 1462 120 0 0 0 0 0 0 51 0 0 0 54 0 116 0 0 15 0 0 179 0 0 105 0 0 0 0 0 152 0 0 0 0 0 0 0 81 0 0 0 0 90 304 116 0 0 0 0 0 0 0 0 112 0 879 0 29 0 0 138 0 0 0 0 129 0 0 0 0 88 0 0 0 0 0 0 0 0 0 196 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 58 0 97 0 0 0 0 0 0 0 0 65 0 0 0 0 0 62 0 93 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 82 0 170 0 0 0 0 0 0 58 0 0 0 84 62 77 0 0 117 138 0 0 0 0 0 0 0 0 0 0 0 0 0 0 32 0 0 0 34 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 60 0 96 0 0 0 0 0 0 0 0 0 0 0 129 0 0 57 0 0 24 0 0 90 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 97 0 0 124 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 86 0 0 0 0 0 0 0 0 0 0 0 0 0 0 22 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 57 0 0 0 0 0 52 151 0 0 0 0 0 0 0 0 0 0 8 0 0 0 0 0 75 0 0 0 0 0 0 0 55 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 24 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
HCB03-V2 0 0 664 0 92 50 559 0 0 63 0 678 385 481 53 165 0 0 326 138 0 73 320 170 0 0 0 0 0 364 305 187 0 0 0 0 311 0 0 0 0 223 0 0 0 0 115 0 42 0 291 139 0 0 133 0 126 0 0 257 0 2237 0 156 0 0 0 0 0 0 0 0 0 227 269 111 2098 230 0 0 123 0 1841 220 170 0 0 0 95 0 113 0 91 0 0 268 177 0 83 72 67 0 481 0 179 0 1772 135 0 250 74 0 0 0 93 0 211 0 0 156 0 0 0 659 554 0 97 0 315 0 107 0 46 0 0 151 0 0 113 1967 0 0 169 0 624 0 0 89 0 0 434 0 0 59 0 0 0 163 108 67 0 0 0 77 0 161 0 0 0 0 0 0 0 0 76 82 0 339 0 0 0 0 0 0 487 140 0 0 0 0 324 155 0 0 0 0 0 96 0 0 0 0 0 0 0 0 0 0 0 0 0 276 0 0 0 0 0 0 0 0 0 0 135 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 106 37 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 104 312 0 0 0 43 0 0 0 280 0 0 0 0 0 72 0 0 0 0 91 0 0 106 24 0 0 0 0 0 166 0 0 0 0 0 0 0 0 0 0 80 0 0 0 0 0 0 0 0 39 0 0 0 130 0 0 0 180 0 92 0 0 0 0 0 64 0 63 62 0 51 0 0 0 0 0 0 0 0 0 0 95 0 0 0 0 144 0 0 0 0 0 0 0 0 0 0 0 0 0 0 154 0 13 0 0 65 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 86 0 117 373 0 0 0 140 0 0 0 0 134 0 0 0 88 0 0 0 0 0 51 0 0 0 0 0 0 0 0 0 0 84 0 0 0 0 0 0 123 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 83 0 0 0 0 0 0 0 0 0 52 0 0 0 0 0 0 136 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 88 0 0 0 0 0 0 8 0 0 0 0 26 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 4 0 0 10 0
HCB04-V2 2261 2120 3782 1573 465 503 3665 0 1219 342 0 2588 0 2120 285 702 0 887 0 636 0 215 0 502 174 198 212 0 0 0 0 414 564 0 144 563 0 225 0 0 807 0 0 180 177 0 319 0 89 206 0 0 0 679 282 0 0 84 149 1018 0 0 0 0 118 0 154 0 141 46 0 0 108 0 0 0 0 0 0 0 0 0 0 0 184 0 0 0 0 242 535 144 0 107 0 0 0 0 0 0 0 0 0 229 0 116 0 157 0 0 439 0 0 0 0 0 0 0 0 0 0 174 0 0 0 0 0 0 0 0 0 0 94 82 0 0 0 0 0 0 0 0 0 0 0 0 19 173 0 0 0 0 0 325 68 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 306 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 74 0 0 0 0 0 0 0 0 172 121 0 0 0 0 0 0 0 0 0 0 0 0 0 0 6 0 0 0 0 0 0 0 0 100 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 73 0 0 0 0 0 0 0 0 0 0 0 0 0 83 0 0 0 0 26 0 0 0 0 0 0 89 0 0 0 255 0 0 0 0 0 0 0 58 0 0 0 0 0 0 0 0 0 242 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 112 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 146 0 0 62 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 54 106 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 113 0 0 0 0 0 51 0 0 0 0 0 0 0 89 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 47 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 15 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 6 0 0 0 13 0 0 0 6 0 0 0 0

#taxonomy table

head(tax_table(psd5)) %>%
  kable(format = "html", col.names = colnames(tax_table(psd5))) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "320px")
Kingdom Phylum Class Order Family Genus
ASV1 Bacteria Firmicutes Bacilli Lactobacillales Streptococcaceae Streptococcus
ASV2 Bacteria Firmicutes Bacilli Lactobacillales Streptococcaceae Streptococcus
ASV3 Bacteria Firmicutes Clostridia Lachnospirales Lachnospiraceae Blautia
ASV4 Bacteria Firmicutes Bacilli Lactobacillales Streptococcaceae Streptococcus
ASV5 Bacteria Proteobacteria Gammaproteobacteria Enterobacterales Enterobacteriaceae Escherichia-Shigella
ASV6 Bacteria Proteobacteria Gammaproteobacteria Enterobacterales Enterobacteriaceae Escherichia-Shigella

#metadata table

as(sample_data(psd5), "data.frame") -> metad
 
metad %>%
  kable(format = "html", col.names = colnames(metad)) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "400px")
SampleID Edad_a_la_inclusion Sexo Respuesta_100 aloTPH_previo_si_no Citogenética_alto_riesgo Afect_Extramedular_PET_screening_si_no Dias_producción CRS_si_no Grado_máx_CRS HOSPITAL Visit
CUN02-V2 CUN02-V2 43 MUJER sCR NO NO SI 120 SI 10 CUN V2
CUN03-V2 CUN03-V2 45 HOMBRE sCR NO SI NO 100 SI 20 CUN V2
HCB01-V2 HCB01-V2 66 HOMBRE PR NO NO SI 110 SI 20 HCB V2
HCB02-V2 HCB02-V2 66 HOMBRE PR NO SI NO 110 SI 20 HCB V2
HCB03-V2 HCB03-V2 60 HOMBRE sCR NO SI SI 110 SI 10 HCB V2
HCB04-V2 HCB04-V2 55 MUJER sCR NO NO NO 110 SI 10 HCB V2
HCB06-V2 HCB06-V2 73 MUJER sCR NO NO NO 90 SI 10 HCB V2
HCB08-V2 HCB08-V2 62 MUJER PR NO NO SI 100 SI 10 HCB V2
HCB09-V2 HCB09-V2 53 HOMBRE sCR NO NO NO 90 SI 10 HCB V2
HCB10-V2 HCB10-V2 52 HOMBRE sCR SI NO SI 100 SI 10 HCB V2
HCB12-V2 HCB12-V2 64 HOMBRE sCR NO SI NO 110 NO NA HCB V2
HVA01-V2 HVA01-V2 58 HOMBRE PR NO NO NO 120 SI 20 HVA V2
HVA02-V2 HVA02-V2 67 HOMBRE sCR NO NO SI 100 SI 20 HVA V2
HVA03-V2 HVA03-V2 61 HOMBRE sCR SI NO SI 140 SI 10 HVA V2
HVA04-V2 HVA04-V2 74 MUJER sCR NO NO NO 100 SI 10 HVA V2
HVR01-V2 HVR01-V2 65 HOMBRE sCR NO SI 110 SI 10 HVR V2
HVR02-V2 HVR02-V2 46 HOMBRE sCR NO SI 120 SI 10 HVR V2
ISBS01-V2 ISBS01-V2 56 HOMBRE sCR NO SI SI 120 SI 10 ISBS V2
ISBS03-V2 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2
ISBS06-V2 ISBS06-V2 60 HOMBRE PR NO NO NO 130 SI 10 ISBS V2
ISBS08-V2 ISBS08-V2 53 HOMBRE sCR NO SI 110 SI 20 ISBS V2
ISBS09-V2 ISBS09-V2 65 MUJER PR NO NO SI 100 SI 10 ISBS V2
ISBS10-V2 ISBS10-V2 68 HOMBRE sCR NO SI SI 140 NO NA ISBS V2

#phylogenetic tree

plot_tree(psd5, method = "treeonly", ladderize = "left")

##Global summary

library(phyloseq)
library(microbiome)
summarize_phyloseq(psd5)
## Compositional = NO2
## 1] Min. number of reads = 146862] Max. number of reads = 1009443] Total number of reads = 7463814] Average number of reads = 32451.3478260875] Median number of reads = 329137] Sparsity = 0.7430333307314036] Any OTU sum to 1 or less? NO8] Number of singletons = 09] Percent of OTUs that are singletons 
##         (i.e. exactly one read detected across all samples)010] Number of sample variables are: 12SampleIDEdad_a_la_inclusionSexoRespuesta_100aloTPH_previo_si_noCitogenética_alto_riesgoAfect_Extramedular_PET_screening_si_noDias_producciónCRS_si_noGrado_máx_CRSHOSPITALVisit2
## [[1]]
## [1] "1] Min. number of reads = 14686"
## 
## [[2]]
## [1] "2] Max. number of reads = 100944"
## 
## [[3]]
## [1] "3] Total number of reads = 746381"
## 
## [[4]]
## [1] "4] Average number of reads = 32451.347826087"
## 
## [[5]]
## [1] "5] Median number of reads = 32913"
## 
## [[6]]
## [1] "7] Sparsity = 0.743033330731403"
## 
## [[7]]
## [1] "6] Any OTU sum to 1 or less? NO"
## 
## [[8]]
## [1] "8] Number of singletons = 0"
## 
## [[9]]
## [1] "9] Percent of OTUs that are singletons \n        (i.e. exactly one read detected across all samples)0"
## 
## [[10]]
## [1] "10] Number of sample variables are: 12"
## 
## [[11]]
##  [1] "SampleID"                              
##  [2] "Edad_a_la_inclusion"                   
##  [3] "Sexo"                                  
##  [4] "Respuesta_100"                         
##  [5] "aloTPH_previo_si_no"                   
##  [6] "Citogenética_alto_riesgo"              
##  [7] "Afect_Extramedular_PET_screening_si_no"
##  [8] "Dias_producción"                       
##  [9] "CRS_si_no"                             
## [10] "Grado_máx_CRS"                         
## [11] "HOSPITAL"                              
## [12] "Visit"

#group all tables

df <- psmelt(psd5)

head(df) %>%
  kable(format = "html", col.names = colnames(df)) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "400px")
OTU Sample Abundance SampleID Edad_a_la_inclusion Sexo Respuesta_100 aloTPH_previo_si_no Citogenética_alto_riesgo Afect_Extramedular_PET_screening_si_no Dias_producción CRS_si_no Grado_máx_CRS HOSPITAL Visit Kingdom Phylum Class Order Family Genus
4871 ASV29 ISBS03-V2 4933 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2 Bacteria Verrucomicrobiota Verrucomicrobiae Verrucomicrobiales Akkermansiaceae Akkermansia
7426 ASV39 ISBS03-V2 4171 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2 Bacteria Verrucomicrobiota Verrucomicrobiae Verrucomicrobiales Akkermansiaceae Akkermansia
7702 ASV40 ISBS03-V2 4020 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2 Bacteria Verrucomicrobiota Verrucomicrobiae Verrucomicrobiales Akkermansiaceae Akkermansia
9466 ASV47 ISBS03-V2 3935 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2 Bacteria Actinobacteriota Actinobacteria Bifidobacteriales Bifidobacteriaceae Bifidobacterium
5117 ASV3 HCB04-V2 3782 HCB04-V2 55 MUJER sCR NO NO NO 110 SI 10 HCB V2 Bacteria Firmicutes Clostridia Lachnospirales Lachnospiraceae Blautia
12070 ASV7 HCB04-V2 3665 HCB04-V2 55 MUJER sCR NO NO NO 110 SI 10 HCB V2 Bacteria Firmicutes Clostridia Lachnospirales Lachnospiraceae Blautia

##Distribution of the samples according to the metadata <> and <>

res <- plot_frequencies(sample_data(psd5), "Respuesta_100", "Sexo")
print(res$plot)
## <environment: 0x55d58a33a180>
# En formato tabla
res$data %>%
  kable(format = "html", col.names = colnames(res$data), digits = 2) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "300px")
Groups Factor n pct .group
PR HOMBRE 4 66.67 1
PR MUJER 2 33.33 1
sCR HOMBRE 12 70.59 2
sCR MUJER 5 29.41 2

#Taxa

taxa_sums(psd5)
##   ASV1   ASV2   ASV3   ASV4   ASV5   ASV6   ASV7   ASV8   ASV9  ASV10  ASV11 
##  16001  14310  12410  10325   8156   7768  10920   6885   7832   6216   6528 
##  ASV12  ASV13  ASV14  ASV15  ASV16  ASV17  ASV18  ASV19  ASV20  ASV21  ASV22 
##   8577   8130   7045   4841   7658   4986   5973   7514   6463   3963   4877 
##  ASV23  ASV24  ASV25  ASV26  ASV27  ASV28  ASV29  ASV30  ASV31  ASV32  ASV33 
##   5876   5391   4502   5060   3177   4933   6493   3813   4751   4545   3698 
##  ASV34  ASV35  ASV36  ASV37  ASV38  ASV39  ASV40  ASV41  ASV42  ASV43  ASV44 
##   4569   3529   3521   3332   5288   5702   5637   1072   3452   2834   3858 
##  ASV45  ASV46  ASV47  ASV48  ASV49  ASV50  ASV51  ASV52  ASV53  ASV54  ASV55 
##   4591   5149   4864   3673   2666   3375   2728   2836   3312    790   4453 
##  ASV56  ASV57  ASV58  ASV59  ASV60  ASV61  ASV62  ASV63  ASV64  ASV65  ASV66 
##   2917   2955   3666   3608   2976   3338   3337    895   3947   1999   2945 
##  ASV67  ASV68  ASV69  ASV70  ASV71  ASV72  ASV73  ASV74  ASV75  ASV76  ASV77 
##   3708    978   2811   2929   2802   2219   2466   2719   2789   2544   2847 
##  ASV78  ASV79  ASV80  ASV81  ASV82  ASV83  ASV84  ASV85  ASV86  ASV87  ASV88 
##   2556   1805   2649   2276   2034   2762   2174   3522   2308   2454   1718 
##  ASV89  ASV90  ASV91  ASV92  ASV93  ASV94  ASV95  ASV96  ASV97  ASV98  ASV99 
##   2236   1975   2689   2876   2240   2675   2609   2291   2035   2622   3095 
## ASV100 ASV101 ASV102 ASV103 ASV104 ASV105 ASV106 ASV107 ASV108 ASV109 ASV110 
##   2188   2879   2822   2208   1884   2089   2419   2361   2883   2181   2027 
## ASV111 ASV112 ASV113 ASV114 ASV115 ASV116 ASV117 ASV118 ASV119 ASV120 ASV121 
##   2262   2030   1472   2159   1794    604   1872   2710   1977   1938   2244 
## ASV122 ASV123 ASV124 ASV125 ASV126 ASV127 ASV128 ASV129 ASV130 ASV131 ASV132 
##   1437   1796   1399   1424   1923   1731   1890   1624   2010   2427   1658 
## ASV133 ASV134 ASV135 ASV136 ASV137 ASV138 ASV139 ASV140 ASV141 ASV142 ASV143 
##   1685   1927   1616   1365   2071   1738   1619   2087   1682   1331   1698 
## ASV144 ASV145 ASV146 ASV147 ASV148 ASV149 ASV150 ASV151 ASV152 ASV153 ASV154 
##   1830   1299   2257   1729   1157   1938   1722   1478   2136   1593   1551 
## ASV155 ASV156 ASV157 ASV158 ASV159 ASV160 ASV161 ASV162 ASV163 ASV164 ASV165 
##   1506   1932   1472   1405   1950   1312   1842   2161    245   1310   1885 
## ASV166 ASV167 ASV168 ASV169 ASV170 ASV171 ASV172 ASV173 ASV174 ASV175 ASV176 
##   1255   1699   1639    781    710   1575   1572   1849   1417   1260   1379 
## ASV177 ASV178 ASV179 ASV180 ASV181 ASV182 ASV183 ASV184 ASV185 ASV186 ASV187 
##   1244   1303   1896   1577   1249   1758    958   1137    951   1170   1079 
## ASV188 ASV189 ASV190 ASV191 ASV192 ASV193 ASV194 ASV195 ASV196 ASV197 ASV198 
##   1294   1326   1470   1223   1042   1733   1057   1014   1705   1231   1483 
## ASV199 ASV200 ASV201 ASV202 ASV203 ASV204 ASV205 ASV206 ASV207 ASV208 ASV209 
##    167   1158   1308   1371   1176   1113   1299    629   1045    628    678 
## ASV210 ASV211 ASV212 ASV213 ASV214 ASV215 ASV216 ASV217 ASV218 ASV219 ASV220 
##    615    697    981   1096    496   1072    923   1114    525    843   1174 
## ASV221 ASV222 ASV223 ASV224 ASV225 ASV226 ASV227 ASV228 ASV229 ASV230 ASV231 
##   1057    569   1058    809    435   1348    775   1176   1314    496    776 
## ASV232 ASV233 ASV234 ASV235 ASV236 ASV237 ASV238 ASV239 ASV240 ASV241 ASV242 
##   1282    901    494    700    769    874    407    712    912    826    562 
## ASV243 ASV244 ASV245 ASV246 ASV247 ASV248 ASV249 ASV250 ASV251 ASV252 ASV253 
##    484   1036    698    506    465   1179    778   1111    438    843   1004 
## ASV254 ASV255 ASV256 ASV257 ASV258 ASV259 ASV260 ASV261 ASV262 ASV263 ASV264 
##    645    760    636   1035    742    861    699    816    856    769    749 
## ASV265 ASV266 ASV267 ASV268 ASV269 ASV270 ASV271 ASV272 ASV273 ASV274 ASV275 
##    843    383    820    687    475    458    548   1043    569    611    499 
## ASV276 ASV277 ASV278 ASV279 ASV280 ASV281 ASV282 ASV283 ASV284 ASV285 ASV286 
##    436    648    988    919    728    365    604    448    547    654    938 
## ASV287 ASV288 ASV289 ASV290 ASV291 ASV292 ASV293 ASV294 ASV295 ASV296 ASV297 
##    542    579    651    526    418    360    538    382    913    326    709 
## ASV298 ASV299 ASV300 ASV301 ASV302 ASV303 ASV304 ASV305 ASV306 ASV307 ASV308 
##    264    355    656    686    544    582    753    501    689    579    406 
## ASV309 ASV310 ASV311 ASV312 ASV313 ASV314 ASV315 ASV316 ASV317 ASV318 ASV319 
##    492    354    725    372    661    317    574    303    470    683    682 
## ASV320 ASV321 ASV322 ASV323 ASV324 ASV325 ASV326 ASV327 ASV328 ASV329 ASV330 
##    650    423    540    719    453    408    603    471    422    142    686 
## ASV331 ASV332 ASV333 ASV334 ASV335 ASV336 ASV337 ASV338 ASV339 ASV340 ASV341 
##    225    305    456    647    499    424    270    247    599    424    466 
## ASV342 ASV343 ASV344 ASV345 ASV346 ASV347 ASV348 ASV349 ASV350 ASV351 ASV352 
##    445    463    325    562    585    334    335    295    559    426    660 
## ASV353 ASV354 ASV355 ASV356 ASV357 ASV358 ASV359 ASV360 ASV361 ASV362 ASV363 
##    507    536    562    583    144    633    455    540    425    223    462 
## ASV364 ASV365 ASV366 ASV367 ASV368 ASV369 ASV370 ASV371 ASV372 ASV373 ASV374 
##    347    501    536    406    405    511    314    415    301    395    341 
## ASV375 ASV376 ASV377 ASV378 ASV379 ASV380 ASV381 ASV382 ASV383 ASV384 ASV385 
##    158    269    404    226    252    396    543    457    531    385    301 
## ASV386 ASV387 ASV388 ASV389 ASV390 ASV391 ASV392 ASV393 ASV394 ASV395 ASV396 
##    521    329    385    345    503    502    413    493    343    357    488 
## ASV397 ASV398 ASV399 ASV400 ASV401 ASV402 ASV403 ASV404 ASV405 ASV406 ASV407 
##    137    139    296    305    401    394    462    318    306    279    233 
## ASV408 ASV409 ASV410 ASV411 ASV412 ASV413 ASV414 ASV415 ASV416 ASV417 ASV418 
##    330    227    192    444    281    436    350    282    430    340    330 
## ASV419 ASV420 ASV421 ASV422 ASV423 ASV424 ASV425 ASV426 ASV427 ASV428 ASV429 
##    427    297    358    336    256    307    161    229    349    325    193 
## ASV430 ASV431 ASV432 ASV433 ASV434 ASV435 ASV436 ASV437 ASV438 ASV439 ASV440 
##    380    273    344    392    245    284    277    163    356    260    248 
## ASV441 ASV442 ASV443 ASV444 ASV445 ASV446 ASV447 ASV448 ASV449 ASV450 ASV451 
##    369    259    187    284    279    359    359    250    356    355    126 
## ASV452 ASV453 ASV454 ASV455 ASV456 ASV457 ASV458 ASV459 ASV460 ASV461 ASV462 
##    349    179    347    224    208    249    220    215    195    263    215 
## ASV463 ASV464 ASV465 ASV466 ASV467 ASV468 ASV469 ASV470 ASV471 ASV472 ASV473 
##    229    330    245    245    299    197    323    265    264    213    200 
## ASV474 ASV475 ASV476 ASV477 ASV478 ASV479 ASV480 ASV481 ASV482 ASV483 ASV484 
##     85    176    253    182    308    246    300     99    131    294    262 
## ASV485 ASV486 ASV487 ASV488 ASV489 ASV490 ASV491 ASV492 ASV493 ASV494 ASV495 
##    286    284    284    278    246    185    272    159    266     80    263 
## ASV496 ASV497 ASV498 ASV499 ASV500 ASV501 ASV502 ASV503 ASV504 ASV505 ASV506 
##    178    257     90    250    205    143    188    245     91    124    159 
## ASV507 ASV508 ASV509 ASV510 ASV511 ASV512 ASV513 ASV514 ASV515 ASV516 ASV517 
##    236    201    233    169     88    228     71    221    179     18    175 
## ASV518 ASV519 ASV520 ASV521 ASV522 ASV523 ASV524 ASV525 ASV526 ASV527 ASV528 
##    206    129    201    138    178    193    132    189    106    188    188 
## ASV529 ASV530 ASV531 ASV532 ASV533 ASV534 ASV535 ASV536 ASV537 ASV538 ASV539 
##    156    162    147    154    154    147    120     65    128    121     99 
## ASV540 ASV541 ASV542 ASV543 ASV544 ASV545 ASV546 ASV547 ASV548 ASV549 ASV550 
##    108    104    103     68    100     26     90     55     75     80     65 
## ASV551 ASV552 ASV553 ASV554 ASV555 ASV556 ASV557 
##     54     50     42     36     35     34     23
kable(tax_table(psd5)[,6], format = "markdown")
Genus
ASV1 Streptococcus
ASV2 Streptococcus
ASV3 Blautia
ASV4 Streptococcus
ASV5 Escherichia-Shigella
ASV6 Escherichia-Shigella
ASV7 Blautia
ASV8 Ruminococcus
ASV9 Streptococcus
ASV10 Escherichia-Shigella
ASV11 Ruminococcus
ASV12 Blautia
ASV13 Collinsella
ASV14 Blautia
ASV15 Escherichia-Shigella
ASV16 Bifidobacterium
ASV17 Ruminococcus
ASV18 Streptococcus
ASV19 Collinsella
ASV20 Bifidobacterium
ASV21 Ruminococcus
ASV22 Streptococcus
ASV23 Collinsella
ASV24 Bifidobacterium
ASV25 Streptococcus
ASV26 Anaerostipes
ASV27 Escherichia-Shigella
ASV28 Ruminococcus
ASV29 Akkermansia
ASV30 Subdoligranulum
ASV31 Collinsella
ASV32 Bifidobacterium
ASV33 Streptococcus
ASV34 Ruminococcus
ASV35 Streptococcus
ASV36 Streptococcus
ASV37 Subdoligranulum
ASV38 Streptococcus
ASV39 Akkermansia
ASV40 Akkermansia
ASV41 Streptococcus
ASV42 Blautia
ASV43 Ruminococcus
ASV44 Anaerostipes
ASV45 Streptococcus
ASV46 Akkermansia
ASV47 Bifidobacterium
ASV48 Ligilactobacillus
ASV49 Streptococcus
ASV50 Anaerostipes
ASV51 Subdoligranulum
ASV52 Blautia
ASV53 Ruminococcus
ASV54 Streptococcus
ASV55 Bifidobacterium
ASV56 Bacteroides
ASV57 Blautia
ASV58 Streptococcus
ASV59 Streptococcus
ASV60 Blautia
ASV61 Ligilactobacillus
ASV62 Holdemanella
ASV63 HT002
ASV64 Blautia
ASV65 Escherichia-Shigella
ASV66 Bacteroides
ASV67 Streptococcus
ASV68 HT002
ASV69 Anaerostipes
ASV70 Enterococcus
ASV71 Enterococcus
ASV72 Bacteroides
ASV73 [Ruminococcus] gnavus group
ASV74 [Eubacterium] hallii group
ASV75 [Eubacterium] hallii group
ASV76 Dorea
ASV77 Holdemanella
ASV78 Bifidobacterium
ASV79 Ruminococcus
ASV80 Ruminococcus
ASV81 Blautia
ASV82 Bacteroides
ASV83 Holdemanella
ASV84 Subdoligranulum
ASV85 Bifidobacterium
ASV86 Bacteroides
ASV87 Blautia
ASV88 Ruminococcus
ASV89 Dorea
ASV90 Erysipelatoclostridium
ASV91 NA
ASV92 Streptococcus
ASV93 Subdoligranulum
ASV94 [Ruminococcus] gnavus group
ASV95 Ligilactobacillus
ASV96 [Eubacterium] hallii group
ASV97 Blautia
ASV98 Agathobacter
ASV99 Blautia
ASV100 Subdoligranulum
ASV101 Bacteroides
ASV102 Bacteroides
ASV103 NA
ASV104 Erysipelatoclostridium
ASV105 Bifidobacterium
ASV106 Streptococcus
ASV107 Holdemanella
ASV108 Bifidobacterium
ASV109 Enterococcus
ASV110 Bifidobacterium
ASV111 NA
ASV112 Faecalibacterium
ASV113 Bacteroides
ASV114 Agathobacter
ASV115 Subdoligranulum
ASV116 HT002
ASV117 [Eubacterium] hallii group
ASV118 Akkermansia
ASV119 Agathobacter
ASV120 [Eubacterium] hallii group
ASV121 Bacteroides
ASV122 Erysipelatoclostridium
ASV123 Blautia
ASV124 Bifidobacterium
ASV125 Bifidobacterium
ASV126 Ruminococcus
ASV127 Dorea
ASV128 Ligilactobacillus
ASV129 NA
ASV130 Faecalibacterium
ASV131 Blautia
ASV132 Bacteroides
ASV133 Fusicatenibacter
ASV134 Streptococcus
ASV135 Romboutsia
ASV136 Subdoligranulum
ASV137 Dorea
ASV138 Enterococcus
ASV139 Subdoligranulum
ASV140 Dialister
ASV141 Faecalibacterium
ASV142 Bacteroides
ASV143 Bifidobacterium
ASV144 Faecalibacterium
ASV145 Bifidobacterium
ASV146 Akkermansia
ASV147 Streptococcus
ASV148 Erysipelatoclostridium
ASV149 Bacteroides
ASV150 Blautia
ASV151 NA
ASV152 Bacteroides
ASV153 [Ruminococcus] torques group
ASV154 NA
ASV155 [Ruminococcus] gnavus group
ASV156 Blautia
ASV157 Fusicatenibacter
ASV158 [Eubacterium] hallii group
ASV159 Blautia
ASV160 Blautia
ASV161 Bacteroides
ASV162 Bacteroides
ASV163 Lacticaseibacillus
ASV164 Bacteroides
ASV165 Bacteroides
ASV166 [Eubacterium] hallii group
ASV167 Dorea
ASV168 Agathobacter
ASV169 Lactobacillus
ASV170 Lactobacillus
ASV171 Streptococcus
ASV172 Blautia
ASV173 Bacteroides
ASV174 Blautia
ASV175 Dorea
ASV176 NA
ASV177 Fusicatenibacter
ASV178 NA
ASV179 Klebsiella
ASV180 Bacteroides
ASV181 [Ruminococcus] torques group
ASV182 Akkermansia
ASV183 Faecalibacterium
ASV184 Faecalibacterium
ASV185 Bifidobacterium
ASV186 Blautia
ASV187 Bacteroides
ASV188 Ligilactobacillus
ASV189 Blautia
ASV190 Dorea
ASV191 Coprococcus
ASV192 [Eubacterium] hallii group
ASV193 Klebsiella
ASV194 Romboutsia
ASV195 Faecalibacterium
ASV196 Akkermansia
ASV197 Ruminococcus
ASV198 Bacteroides
ASV199 Lacticaseibacillus
ASV200 [Ruminococcus] gnavus group
ASV201 [Clostridium] innocuum group
ASV202 Bacteroides
ASV203 Faecalibacterium
ASV204 Ruminococcus
ASV205 Streptococcus
ASV206 Alistipes
ASV207 Lachnoclostridium
ASV208 Limosilactobacillus
ASV209 Enterococcus
ASV210 Enterococcus
ASV211 Fusicatenibacter
ASV212 Coprococcus
ASV213 Blautia
ASV214 Alistipes
ASV215 Streptococcus
ASV216 [Clostridium] innocuum group
ASV217 Dorea
ASV218 Limosilactobacillus
ASV219 [Ruminococcus] torques group
ASV220 Anaerostipes
ASV221 Streptococcus
ASV222 NA
ASV223 Blautia
ASV224 Alistipes
ASV225 Lactobacillus
ASV226 Bacteroides
ASV227 Alistipes
ASV228 Bacteroides
ASV229 NA
ASV230 NA
ASV231 Faecalibacterium
ASV232 NA
ASV233 Ligilactobacillus
ASV234 Enterococcus
ASV235 Eggerthella
ASV236 Bacteroides
ASV237 [Ruminococcus] torques group
ASV238 Lactobacillus
ASV239 Romboutsia
ASV240 Blautia
ASV241 Streptococcus
ASV242 Faecalibacterium
ASV243 Collinsella
ASV244 UCG-002
ASV245 [Ruminococcus] torques group
ASV246 Erysipelatoclostridium
ASV247 Limosilactobacillus
ASV248 Klebsiella
ASV249 Alistipes
ASV250 Streptococcus
ASV251 Collinsella
ASV252 Streptococcus
ASV253 Blautia
ASV254 Romboutsia
ASV255 Streptococcus
ASV256 Parabacteroides
ASV257 Methanobrevibacter
ASV258 UBA1819
ASV259 [Ruminococcus] torques group
ASV260 Coprococcus
ASV261 UCG-002
ASV262 Blautia
ASV263 NA
ASV264 Hespellia
ASV265 Anaerostipes
ASV266 Enterococcus
ASV267 Roseburia
ASV268 Coprococcus
ASV269 Ruminococcus
ASV270 NA
ASV271 Streptococcus
ASV272 NA
ASV273 Blautia
ASV274 Monoglobus
ASV275 Sellimonas
ASV276 Faecalibacterium
ASV277 [Eubacterium] hallii group
ASV278 Streptococcus
ASV279 Senegalimassilia
ASV280 Bacteroides
ASV281 Collinsella
ASV282 Bifidobacterium
ASV283 Eggerthella
ASV284 Lachnoclostridium
ASV285 UBA1819
ASV286 Bacteroides
ASV287 Lachnoclostridium
ASV288 Parabacteroides
ASV289 Butyricicoccus
ASV290 Veillonella
ASV291 Erysipelatoclostridium
ASV292 Collinsella
ASV293 Erysipelatoclostridium
ASV294 NA
ASV295 Bifidobacterium
ASV296 Alistipes
ASV297 Bacteroides
ASV298 Limosilactobacillus
ASV299 Faecalibacterium
ASV300 Bifidobacterium
ASV301 Bacteroides
ASV302 [Clostridium] innocuum group
ASV303 Haemophilus
ASV304 Streptococcus
ASV305 Sellimonas
ASV306 Roseburia
ASV307 Alistipes
ASV308 Eggerthella
ASV309 Monoglobus
ASV310 Lachnoclostridium
ASV311 Bacteroides
ASV312 Lachnoclostridium
ASV313 Dorea
ASV314 NA
ASV315 Streptococcus
ASV316 NA
ASV317 NA
ASV318 UCG-002
ASV319 Hespellia
ASV320 Bacteroides
ASV321 [Eubacterium] hallii group
ASV322 Clostridium sensu stricto 1
ASV323 Methanobrevibacter
ASV324 Alistipes
ASV325 Anaerostipes
ASV326 Bacteroides
ASV327 Monoglobus
ASV328 Parabacteroides
ASV329 Subdoligranulum
ASV330 Senegalimassilia
ASV331 Erysipelatoclostridium
ASV332 NA
ASV333 Alistipes
ASV334 Coprococcus
ASV335 [Eubacterium] hallii group
ASV336 Streptococcus
ASV337 Veillonella
ASV338 NA
ASV339 Streptococcus
ASV340 Coprococcus
ASV341 Bifidobacterium
ASV342 Lachnoclostridium
ASV343 [Clostridium] innocuum group
ASV344 Lactococcus
ASV345 Bacteroides
ASV346 Blautia
ASV347 Bacteroides
ASV348 Lachnoclostridium
ASV349 NA
ASV350 UCG-002
ASV351 Veillonella
ASV352 Haemophilus
ASV353 Roseburia
ASV354 Coprococcus
ASV355 Methanobrevibacter
ASV356 Streptococcus
ASV357 Lactococcus
ASV358 Bacteroides
ASV359 Christensenellaceae R-7 group
ASV360 Roseburia
ASV361 Anaerostipes
ASV362 NA
ASV363 UBA1819
ASV364 Faecalibacterium
ASV365 Bacteroides
ASV366 Senegalimassilia
ASV367 [Ruminococcus] torques group
ASV368 [Eubacterium] hallii group
ASV369 Blautia
ASV370 Veillonella
ASV371 Bacteroides
ASV372 [Eubacterium] hallii group
ASV373 Bacteroides
ASV374 Alistipes
ASV375 Lachnoclostridium
ASV376 Family XIII AD3011 group
ASV377 Bacteroides
ASV378 [Eubacterium] hallii group
ASV379 Eggerthella
ASV380 Subdoligranulum
ASV381 Blautia
ASV382 Bacteroides
ASV383 Coprococcus
ASV384 NA
ASV385 Bifidobacterium
ASV386 [Ruminococcus] torques group
ASV387 Erysipelatoclostridium
ASV388 UCG-005
ASV389 Flavonifractor
ASV390 Bacteroides
ASV391 Bacteroides
ASV392 Subdoligranulum
ASV393 [Ruminococcus] gauvreauii group
ASV394 Christensenellaceae R-7 group
ASV395 NA
ASV396 Streptococcus
ASV397 Lachnoclostridium
ASV398 Sellimonas
ASV399 Anaerostipes
ASV400 Blautia
ASV401 Incertae Sedis
ASV402 Coprococcus
ASV403 Lachnospiraceae NK4A136 group
ASV404 Sellimonas
ASV405 Monoglobus
ASV406 Phascolarctobacterium
ASV407 Incertae Sedis
ASV408 Slackia
ASV409 [Eubacterium] brachy group
ASV410 Bacteroides
ASV411 Christensenellaceae R-7 group
ASV412 Agathobacter
ASV413 Blautia
ASV414 Senegalimassilia
ASV415 NA
ASV416 Hespellia
ASV417 Alistipes
ASV418 Streptococcus
ASV419 Eubacterium
ASV420 UCG-002
ASV421 NA
ASV422 Clostridium sensu stricto 1
ASV423 Collinsella
ASV424 Streptococcus
ASV425 NA
ASV426 Family XIII AD3011 group
ASV427 Bacteroides
ASV428 Butyricicoccus
ASV429 Erysipelatoclostridium
ASV430 Coprococcus
ASV431 NA
ASV432 Veillonella
ASV433 NA
ASV434 Family XIII AD3011 group
ASV435 Bacteroides
ASV436 UCG-005
ASV437 Incertae Sedis
ASV438 Streptococcus
ASV439 Blautia
ASV440 Bacteroides
ASV441 Incertae Sedis
ASV442 [Ruminococcus] gauvreauii group
ASV443 Eggerthella
ASV444 Bacteroides
ASV445 Clostridium sensu stricto 1
ASV446 Coprococcus
ASV447 UBA1819
ASV448 Bacteroides
ASV449 Blautia
ASV450 Blautia
ASV451 [Eubacterium] siraeum group
ASV452 Bacteroides
ASV453 Alistipes
ASV454 Christensenellaceae R-7 group
ASV455 Gemella
ASV456 Intestinibacter
ASV457 Incertae Sedis
ASV458 NA
ASV459 Sellimonas
ASV460 Streptococcus
ASV461 Roseburia
ASV462 Subdoligranulum
ASV463 Bifidobacterium
ASV464 Lachnoclostridium
ASV465 Parabacteroides
ASV466 Blautia
ASV467 Incertae Sedis
ASV468 Roseburia
ASV469 Alistipes
ASV470 Anaerostipes
ASV471 Slackia
ASV472 Christensenellaceae R-7 group
ASV473 NK4A214 group
ASV474 Coprobacillus
ASV475 NA
ASV476 Oscillibacter
ASV477 Faecalibacterium
ASV478 [Eubacterium] eligens group
ASV479 Streptococcus
ASV480 Intestinimonas
ASV481 Latilactobacillus
ASV482 NA
ASV483 Subdoligranulum
ASV484 Blautia
ASV485 Parabacteroides
ASV486 Christensenellaceae R-7 group
ASV487 NA
ASV488 Olsenella
ASV489 Oscillibacter
ASV490 NA
ASV491 Lachnoclostridium
ASV492 Roseburia
ASV493 Parabacteroides
ASV494 Latilactobacillus
ASV495 Family XIII AD3011 group
ASV496 NA
ASV497 Blautia
ASV498 Flavonifractor
ASV499 Bacteroides
ASV500 DTU089
ASV501 Odoribacter
ASV502 CAG-56
ASV503 Christensenellaceae R-7 group
ASV504 Blautia
ASV505 Streptococcus
ASV506 Roseburia
ASV507 UCG-002
ASV508 Incertae Sedis
ASV509 Terrisporobacter
ASV510 Oscillibacter
ASV511 Lachnoclostridium
ASV512 Lachnospiraceae NK4A136 group
ASV513 Lachnoclostridium
ASV514 Anaerostipes
ASV515 Incertae Sedis
ASV516 Solobacterium
ASV517 Alistipes
ASV518 [Ruminococcus] torques group
ASV519 Intestinimonas
ASV520 Christensenellaceae R-7 group
ASV521 [Eubacterium] siraeum group
ASV522 Family XIII AD3011 group
ASV523 Christensenellaceae R-7 group
ASV524 Streptococcus
ASV525 Bacteroides
ASV526 Family XIII AD3011 group
ASV527 Coprobacillus
ASV528 UBA1819
ASV529 Oscillibacter
ASV530 Christensenellaceae R-7 group
ASV531 NK4A214 group
ASV532 Incertae Sedis
ASV533 [Eubacterium] nodatum group
ASV534 UC5-1-2E3
ASV535 Christensenellaceae R-7 group
ASV536 Family XIII UCG-001
ASV537 NA
ASV538 Alistipes
ASV539 [Eubacterium] brachy group
ASV540 Lachnospira
ASV541 Holdemania
ASV542 Bilophila
ASV543 Scardovia
ASV544 Family XIII AD3011 group
ASV545 TM7x
ASV546 NA
ASV547 Anaerotruncus
ASV548 Actinomyces
ASV549 Solobacterium
ASV550 NA
ASV551 Lachnospiraceae FCS020 group
ASV552 Turicibacter
ASV553 Bifidobacterium
ASV554 NA
ASV555 UCG-005
ASV556 Parabacteroides
ASV557 Candidatus Soleaferrea
#ASVs of the same gender at the sequence level are different even though they correspond to the same gender. Visualization of sequences of the taxonomic range in phylogenetic trees

# cluster by gender
psd5.genus = tax_glom(psd5, "Genus", NArm = FALSE)
# Luego por altura en el árbol filogenético
h1 = 0.4
psd5.tip = tip_glom(psd5, h = h1)
# plot a comparison to visualize the differences
multiPlotTitleTextSize = 15
p2tree = plot_tree(psd5, method = "treeonly",
                   ladderize = "left",
                   title = "Sin aglomeración") +
  theme(plot.title = element_text(size = multiPlotTitleTextSize))
p3tree = plot_tree(psd5.genus, method = "treeonly",
                   ladderize = "left", title = "A nivel de género") +
  theme(plot.title = element_text(size = multiPlotTitleTextSize))
p4tree = plot_tree(psd5.tip, method = "treeonly",
                   ladderize = "left", title = "Por altura") +
  theme(plot.title = element_text(size = multiPlotTitleTextSize))

# plot all trees 
grid.arrange(nrow = 1, p2tree, p3tree, p4tree)

####DIVERSITY ANALYSIS####

#Richness Summary Graphics
plot_richness(psd5)
## Warning in estimate_richness(physeq, split = TRUE, measures = measures): The data you have provided does not have
## any singletons. This is highly suspicious. Results of richness
## estimates (for example) are probably unreliable, or wrong, if you have already
## trimmed low-abundance taxa from the data.
## 
## We recommended that you find the un-trimmed data and retry.

#Richness Summary Graphics group by <<Respuesta_100>>
(p = plot_richness(psd5, x = "Respuesta_100"))
## Warning in estimate_richness(physeq, split = TRUE, measures = measures): The data you have provided does not have
## any singletons. This is highly suspicious. Results of richness
## estimates (for example) are probably unreliable, or wrong, if you have already
## trimmed low-abundance taxa from the data.
## 
## We recommended that you find the un-trimmed data and retry.

#Richness Summary Graphics group by <<Respuesta_100>>
plot_richness(psd5, color = "Respuesta_100", x = "Respuesta_100", measures = c("Observed", "Chao1", "Shannon")) + geom_boxplot(aes(fill = Respuesta_100), alpha=.7) + scale_color_manual(values = c("#a6cee3", "#b2df8a")) + scale_fill_manual(values = c("#a6cee3", "#b2df8a"))
## Warning in estimate_richness(physeq, split = TRUE, measures = measures): The data you have provided does not have
## any singletons. This is highly suspicious. Results of richness
## estimates (for example) are probably unreliable, or wrong, if you have already
## trimmed low-abundance taxa from the data.
## 
## We recommended that you find the un-trimmed data and retry.

#Check if there is any significant effect of alpha diversity based on response type

# save a dataframe with the alpha diversity measures
alpha_pd <- estimate_pd(psd5)
## Calculating Faiths PD-index...
# combine the metadata with alpha.diversity
data <- cbind(sample_data(psd5), alpha_pd) 

# calculate an ANOVA
psd5.anova <- aov(PD ~ Respuesta_100, data) 
# install.packages("xtable")
library(xtable)
psd5.anova.table <- xtable(psd5.anova)
kable(psd5.anova.table, caption = "Tabla ANOVA", digits = 5, format = "markdown")
Tabla ANOVA
Df Sum Sq Mean Sq F value Pr(>F)
Respuesta_100 1 1.38598 1.38598 0.0554 0.81619
Residuals 21 525.32500 25.01548 NA NA
#Mann-Whitney-Wilcoxon Test
psd5.wilcox <- wilcox.test(PD ~ Respuesta_100, data) 
psd5.wilcox 
## 
##  Wilcoxon rank sum exact test
## 
## data:  PD by Respuesta_100
## W = 57, p-value = 0.7077
## alternative hypothesis: true location shift is not equal to 0
#Tablle of divsersity scores
tab <- global(psd5, index = "all")
## Warning: The microbiome::global function is deprecated. 
##         Use the function microbiome::alpha instead.
## Observed richness
## Other forms of richness
## Diversity
## Evenness
## Dominance
## Rarity
head(tab) %>%
  kable(format = "html", col.names = colnames(tab), digits = 2) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "300px")
observed chao1 diversity_inverse_simpson diversity_gini_simpson diversity_shannon diversity_fisher diversity_coverage evenness_camargo evenness_pielou evenness_simpson evenness_evar evenness_bulla dominance_dbp dominance_dmn dominance_absolute dominance_relative dominance_simpson dominance_core_abundance dominance_gini rarity_log_modulo_skewness rarity_low_abundance rarity_rare_abundance
CUN02-V2 121 121 44.17 0.98 4.23 15.62 17 0.50 0.88 0.37 0.45 0.57 0.06 0.11 2081 0.06 0.02 0.07 0.90 1.59 0.03 0.24
CUN03-V2 188 188 76.19 0.99 4.82 27.05 37 0.50 0.92 0.41 0.59 0.68 0.05 0.10 1391 0.05 0.01 0.24 0.82 1.74 0.04 0.25
HCB01-V2 193 193 89.62 0.99 4.87 26.84 37 0.51 0.93 0.46 0.57 0.66 0.03 0.07 1192 0.03 0.01 0.16 0.81 1.65 0.06 0.31
HCB02-V2 88 88 29.97 0.97 3.94 12.28 13 0.77 0.88 0.34 0.55 0.61 0.10 0.19 1639 0.10 0.03 0.04 0.92 1.76 0.01 0.58
HCB03-V2 126 126 36.74 0.97 4.20 16.75 15 0.49 0.87 0.29 0.49 0.58 0.07 0.14 2237 0.07 0.03 0.21 0.90 1.67 0.02 0.19
HCB04-V2 83 83 22.14 0.95 3.61 10.14 8 0.33 0.82 0.27 0.32 0.49 0.10 0.20 3782 0.10 0.05 0.66 0.95 1.59 0.01 0.10
library(microbiome)
# Richness
tab <- richness(psd5)
# Dominance
tab <- dominance(psd5, index = "all")
# Rarity
tab <- rarity(psd5, index = "all")
# Coverage
tab <- coverage(psd5, threshold = 0.5)
# inequality
tab <- inequality(psd5)
# evennes
tab <- evenness(psd5, "all")
library(ggpubr)
# Generate a `phyloseq` object without taxa that adds 0 reads

psd5.2 <- prune_taxa(taxa_sums(psd5) > 0, psd5)
# calculate the diversity indices
tab <- diversities(psd5.2, index = "all")
## Warning: 'diversities' is deprecated.
## Use 'diversity' instead.
## See help("Deprecated") and help("The microbiome::diversities function has been 
##         replaced by function microbiome::alpha. Update your code accordingly.-deprecated").
## Observed richness
## Other forms of richness
## Diversity
## Evenness
## Dominance
## Rarity
# visualize the results table
head(tab) %>%
  kable(format = "html", col.names = colnames(tab), digits = 2) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "310px")
observed chao1 diversity_inverse_simpson diversity_gini_simpson diversity_shannon diversity_fisher diversity_coverage evenness_camargo evenness_pielou evenness_simpson evenness_evar evenness_bulla dominance_dbp dominance_dmn dominance_absolute dominance_relative dominance_simpson dominance_core_abundance dominance_gini rarity_log_modulo_skewness rarity_low_abundance rarity_rare_abundance
CUN02-V2 121 121 44.17 0.98 4.23 15.62 17 0.50 0.88 0.37 0.45 0.57 0.06 0.11 2081 0.06 0.02 0.07 0.90 1.59 0.03 0.24
CUN03-V2 188 188 76.19 0.99 4.82 27.05 37 0.50 0.92 0.41 0.59 0.68 0.05 0.10 1391 0.05 0.01 0.24 0.82 1.74 0.04 0.25
HCB01-V2 193 193 89.62 0.99 4.87 26.84 37 0.51 0.93 0.46 0.57 0.66 0.03 0.07 1192 0.03 0.01 0.16 0.81 1.65 0.06 0.31
HCB02-V2 88 88 29.97 0.97 3.94 12.28 13 0.77 0.88 0.34 0.55 0.61 0.10 0.19 1639 0.10 0.03 0.04 0.92 1.76 0.01 0.58
HCB03-V2 126 126 36.74 0.97 4.20 16.75 15 0.49 0.87 0.29 0.49 0.58 0.07 0.14 2237 0.07 0.03 0.21 0.90 1.67 0.02 0.19
HCB04-V2 83 83 22.14 0.95 3.61 10.14 8 0.33 0.82 0.27 0.32 0.49 0.10 0.20 3782 0.10 0.05 0.66 0.95 1.59 0.01 0.10
#extract metadata
psd5.2.meta <- meta(psd5.2)

head(psd5.2.meta) %>%
  kable(format = "html", col.names = colnames(psd5.2.meta), digits = 2) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "500px")
SampleID Edad_a_la_inclusion Sexo Respuesta_100 aloTPH_previo_si_no Citogenética_alto_riesgo Afect_Extramedular_PET_screening_si_no Dias_producción CRS_si_no Grado_máx_CRS HOSPITAL Visit
CUN02-V2 CUN02-V2 43 MUJER sCR NO NO SI 120 SI 10 CUN V2
CUN03-V2 CUN03-V2 45 HOMBRE sCR NO SI NO 100 SI 20 CUN V2
HCB01-V2 HCB01-V2 66 HOMBRE PR NO NO SI 110 SI 20 HCB V2
HCB02-V2 HCB02-V2 66 HOMBRE PR NO SI NO 110 SI 20 HCB V2
HCB03-V2 HCB03-V2 60 HOMBRE sCR NO SI SI 110 SI 10 HCB V2
HCB04-V2 HCB04-V2 55 MUJER sCR NO NO NO 110 SI 10 HCB V2
#add the diversity table to the metadata
psd5.2.meta$Shannon <- tab$diversity_shannon
# get the variables from our `phyloseq` object
rpp <- levels(psd5.2.meta$Respuesta_100)
#Create a list of what we want to compare

pares.rpp <- combn(seq_along(rpp), 2, simplify = FALSE, FUN = function(i)rpp[i])
# Print results
print(pares.rpp)
## [[1]]
## [1] "PR"  "sCR"
#Plot 
p1 <- ggviolin(psd5.2.meta, x = "Respuesta_100", y = "Shannon",
 add = "boxplot", fill = "Respuesta_100", palette = c("#a6cee3", "#b2df8a"))  
print(p1)

#Plot with compare means
p1 <- p1 + stat_compare_means(comparisons = pares.rpp)
print(p1)

###Beta diversity

#UniFrac
psd5.mds.unifrac <- ordinate(psd5, method = "MDS", distance = "unifrac")
## Warning in UniFrac(physeq, ...): Randomly assigning root as -- ASV415 -- in the
## phylogenetic tree in the data you provided.
evals <- psd5.mds.unifrac$values$Eigenvalues
pord1 <- plot_ordination(psd5, psd5.mds.unifrac, color = "Respuesta_100") +
  labs(col = "Respuesta_100") +
  coord_fixed(sqrt(evals[2] / evals[1]))

#Bray-Curtis
psd5.mds.bray <- ordinate(psd5, method = "MDS", distance = "bray")
evals <- psd5.mds.bray$values$Eigenvalues
pord2 <- plot_ordination(psd5, psd5.mds.bray, color = "Respuesta_100") +
  labs(col = "Respuesta_100") +
  coord_fixed(sqrt(evals[2] / evals[1]))

grid.arrange(pord1, pord2)

#DPCoA
psd5.dpcoa.unifrac <- ordinate(psd5, method = "DPCoA", distance = "dpcoa")
evals <- psd5.dpcoa.unifrac$eig
pord3 <- plot_ordination(psd5, psd5.dpcoa.unifrac, color = "Respuesta_100", shape = "Sexo") +
  labs(col = "Respuesta_100") +
  coord_fixed(sqrt(evals[2] / evals[1])) +
  scale_color_manual(values = c("#a6cee3", "#b2df8a")) + 
  scale_fill_manual(values = c("#a6cee3",  "#b2df8a")) +
  geom_point(size=4)
## Warning: non-unique values when setting 'row.names':
## Species coordinates not found directly in ordination object. Attempting weighted average (`vegan::wascores`)
pord3

#multidimensional scaling

library(tsnemicrobiota)
tsne_res <- tsne_phyloseq(psd5, distance= "dpcoa", perplexity = 8, verbose=0, rng_seed = 3901)

# Plot
pord4 <- plot_tsne_phyloseq(psd5, tsne_res, color = "Respuesta_100", shape = "Sexo") +
  geom_point(size=4) +
  scale_color_manual(values = c("#a6cee3",  "#b2df8a")) + 
  scale_fill_manual(values = c("#a6cee3",  "#b2df8a"))

grid.arrange(pord3, pord4)
## Warning: Use of `plot_df[[axes[1]]]` is discouraged. Use `.data[[axes[1]]]`
## instead.
## Warning: Use of `plot_df[[axes[2]]]` is discouraged. Use `.data[[axes[2]]]`
## instead.
## Warning: Use of `plot_df[[axes[1]]]` is discouraged. Use `.data[[axes[1]]]`
## instead.
## Warning: Use of `plot_df[[axes[2]]]` is discouraged. Use `.data[[axes[2]]]`
## instead.

###Abundance analysis

# Necesitamos obtener las taxa más abundantes, en este caso el top 15
top15 <- get_top_taxa(physeq_obj = psd5, n = 15, relative = T,
                       discard_other = T, other_label = "Other")
# Ya que no todas las taxa fueron clasificadas a nivel de especie, generamos etiquetas compuestas de distintos rangos taxonómicos para el gráfico
top15 <- name_taxa(top15, label = "", species = F, other_label = "Other")
# Finalmente graficamos
fantaxtic_bar(top15, color_by = "Family", label_by = "Genus", facet_by = NULL, grid_by = NULL, other_color = "Grey") -> ptop15
##                Level N.color.shades Central.color
## 1   Streptococcaceae              4       #6495ed
## 2    Lachnospiraceae              4       #ff7256
## 3 Enterobacteriaceae              2       #edbc64
## 4  Coriobacteriaceae              3       #8470ff
## 5 Bifidobacteriaceae              2       #8ee5ee
#Barplot of abundance
ptop15

#Groupo by Respuesta_100
fantaxtic_bar(top15, color_by = "Family", label_by = "Genus", facet_by = "Respuesta_100", grid_by = NULL, other_color = "Grey") -> ptop15.2
##                Level N.color.shades Central.color
## 1   Streptococcaceae              4       #6495ed
## 2    Lachnospiraceae              4       #ff7256
## 3 Enterobacteriaceae              2       #edbc64
## 4  Coriobacteriaceae              3       #8470ff
## 5 Bifidobacteriaceae              2       #8ee5ee
#Plot 
ptop15.2

library(ampvis2)
# extract the read counts table and the taxonomy table from the psd5 object
#  generate a copy so as not to overwrite psd5
obj <- psd5
# change the orientation of the otu_table
t(otu_table(obj)) -> otu_table(obj)
# Extraemos las tablas
otutable <- data.frame(OTU = rownames(phyloseq::otu_table(obj)@.Data),
                       phyloseq::otu_table(obj)@.Data,
                       phyloseq::tax_table(obj)@.Data,
                       check.names = FALSE
)
# extract la metadada
metadata <- data.frame(phyloseq::sample_data(obj), 
                       check.names = FALSE
)

#ampvis2 requires that 1) the taxonomic ranges are seven and go from Kingdom to Species and 2) the first column of the metadata is the identifier of each sample
# Then we duplicate the Genus column and rename it to Species
otutable$Species = otutable$Genus
# Reordenamos la metadata
metadata <- metadata[,c("SampleID", "Edad_a_la_inclusion", "Sexo",  "Respuesta_100",    "aloTPH_previo_si_no",  "Citogenética_alto_riesgo", "Afect_Extramedular_PET_screening_si_no",   "Dias_producción",  "CRS_si_no", "Grado_máx_CRS",   "HOSPITAL", "Visit")]

# generate the ampvis object

av2 <- amp_load(otutable, metadata)
#estructura de las comunidades en las curvas de abundancia
amp_rankabundance(av2, log10_x = T, group_by = "Sexo")

#estructura de las comunidades en las curvas de abundancia
amp_rankabundance(av2, log10_x = T, group_by = "Respuesta_100")

#estructura de las comunidades en las curvas de abundancia
amp_rankabundance(av2, log10_x = T, group_by = "HOSPITAL")

#heatmap of abundance
amp_heatmap(av2, 
            group_by = "Respuesta_100", 
            facet_by = "Sexo", 
            plot_values = TRUE,
            tax_show = 20,
            tax_aggregate = "Genus",
            tax_add = "Phylum",
            plot_colorscale = "sqrt",
            plot_legendbreaks = c(1, 5, 10)
            )

#Box Plots.
amp_boxplot(av2,
            group_by = "Respuesta_100",
            tax_show = 20,
            tax_aggregate = "Genus",
            tax_add = "Phylum",
            adjust_zero = T,
            plot_log = T) +
  scale_color_manual(values = c("#a6cee3", "#b2df8a")) + 
  scale_fill_manual(values = c("#a6cee3", "#b2df8a"))

#Venn Diagram
# some of these microorganisms are shared between all samples
amp_venn(av2, group_by = "Respuesta_100", cut_a = 0, cut_f = 50, text_size = 3)

#Differential abundance analysis

#we want to determine exactly which taxa is more represented in one condition versus another and to what extent. DESeq2 package,

# Create a DESeq2 object with the `phyloseq_to_deseq2` function
diagdds = phyloseq_to_deseq2(psd4, ~Respuesta_100)
## converting counts to integer mode
# We calculate the size factors as part of the normalization of the samples
# calculate geometric means prior to estimate size factors
gm_mean = function(x, na.rm=TRUE){
    exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))
}
geoMeans = apply(counts(diagdds), 1, gm_mean)
diagdds = estimateSizeFactors(diagdds, geoMeans = geoMeans)

# We normalize and perform the parametric Wald test to determine differentially abundant taxa.
diagdds = DESeq(diagdds, test="Wald", fitType="local")
## using pre-existing size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
## -- replacing outliers and refitting for 297 genes
## -- DESeq argument 'minReplicatesForReplace' = 7 
## -- original counts are preserved in counts(dds)
## estimating dispersions
## fitting model and testing
# Create a DESeq2 object with the `phyloseq_to_deseq2` function
diagdds = phyloseq_to_deseq2(psd5, ~Respuesta_100)
## converting counts to integer mode
# We calculate the size factors as part of the normalization of the samples
# calculate geometric means prior to estimate size factors
gm_mean = function(x, na.rm=TRUE){
    exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))
}
geoMeans = apply(counts(diagdds), 1, gm_mean)
diagdds = estimateSizeFactors(diagdds, geoMeans = geoMeans)

# We normalize and perform the parametric Wald test to determine differentially abundant taxa.
diagdds = DESeq(diagdds, test="Wald", fitType="local")
## using pre-existing size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
## -- replacing outliers and refitting for 297 genes
## -- DESeq argument 'minReplicatesForReplace' = 7 
## -- original counts are preserved in counts(dds)
## estimating dispersions
## fitting model and testing
## estimating size factors
plotMA(diagdds)

## estimating dispersions
plotDispEsts(diagdds)

# save the results 
res = results(diagdds, cooksCutoff = FALSE)
# reorder 
res = res[order(res$padj, na.last=NA), ]
#especify contrast
res.PR.sCR <- results(diagdds, contrast=c("Respuesta_100","PR","sCR"))
#Visualizate generated table
res.PR.sCR  %>%
  kable(format = "html", col.names = colnames(res.PR.sCR )) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "400px")
baseMean log2FoldChange lfcSE stat pvalue padj
ASV1 729.6119002 -1.4986609 1.481991 -1.0112485 0.3118975 0.9759373
ASV2 647.9485881 -1.7206118 1.752171 -0.9819884 0.3261055 0.9919616
ASV3 567.9817458 -1.6673863 1.140755 -1.4616510 0.1438369 0.8239591
ASV4 469.9120778 -1.8321980 1.676678 -1.0927550 0.2745014 0.9181599
ASV5 335.6251807 0.1718067 1.840638 0.0933409 0.9256328 0.9919616
ASV6 318.8184377 0.2142279 1.845320 0.1160926 0.9075792 0.9919616
ASV7 502.3517870 -1.6004280 1.459003 -1.0969324 0.2726710 0.9181599
ASV8 231.8964389 0.1585241 3.285750 0.0482459 0.9615202 0.9919616
ASV9 354.9656963 -1.5880479 1.624910 -0.9773146 0.3284134 0.9919616
ASV10 256.4548414 0.2538495 1.762246 0.1440489 0.8854619 0.9919616
ASV11 220.4540728 0.2587368 3.285756 0.0787450 0.9372355 0.9919616
ASV12 394.4344611 -1.4672698 1.384096 -1.0600924 0.2891026 0.9506084
ASV13 387.6501461 0.1555353 1.799128 0.0864504 0.9311084 0.9919616
ASV14 324.1069347 -1.1306624 1.612926 -0.7010006 0.4833026 0.9919616
ASV15 200.5535457 0.0647585 1.799807 0.0359808 0.9712977 0.9919616
ASV16 344.9188676 -0.5298251 1.510569 -0.3507453 0.7257794 0.9919616
ASV17 98.1964810 1.9335166 3.286040 0.5884032 0.5562617 0.9919616
ASV18 276.9018480 -1.8157152 1.969764 -0.9217930 0.3566366 0.9919616
ASV19 356.5740209 0.3518504 1.923800 0.1828934 0.8548817 0.9919616
ASV20 287.9218237 -0.3332220 1.845926 -0.1805175 0.8567463 0.9919616
ASV21 79.4269991 2.0999301 3.286195 0.6390157 0.5228127 0.9919616
ASV22 269.4452622 2.7801277 3.056320 0.9096325 0.3630164 0.9919616
ASV23 278.4514632 0.3426924 1.952610 0.1755047 0.8606830 0.9919616
ASV24 243.9912163 -0.4321198 1.982185 -0.2180017 0.8274278 0.9919616
ASV25 203.9225727 4.3251729 3.203453 1.3501595 0.1769648 0.8239591
ASV26 212.9596269 -0.7757047 1.414269 -0.5484845 0.5833593 0.9919616
ASV27 130.9438077 0.2653469 2.032955 0.1305228 0.8961528 0.9919616
ASV28 220.1562885 1.1202239 3.120900 0.3589426 0.7196380 0.9919616
ASV29 77.4021602 2.0629353 3.286211 0.6277551 0.5301644 0.9919616
ASV30 160.9733415 0.9049913 2.271571 0.3983988 0.6903363 0.9919616
ASV31 224.1176959 0.2021397 1.976725 0.1022599 0.9185504 0.9919616
ASV32 205.1326821 -0.3727775 1.886760 -0.1975755 0.8433772 0.9919616
ASV33 171.2487384 -1.9942607 2.121064 -0.9402172 0.3471061 0.9919616
ASV34 203.2785965 0.8817885 3.285741 0.2683682 0.7884159 0.9919616
ASV35 170.9246972 3.9919722 2.560842 1.5588516 NA NA
ASV36 161.5216435 -2.0737110 2.122820 -0.9768663 0.3286354 0.9919616
ASV37 142.5629474 0.7749435 2.295679 0.3375662 0.7356901 0.9919616
ASV38 224.9469085 -1.7253151 2.099813 -0.8216518 0.4112751 0.9919616
ASV39 76.5044086 1.7966353 3.175117 0.5658484 0.5714968 0.9919616
ASV40 80.0219857 2.1775840 3.286198 0.6626453 0.5075577 0.9919616
ASV41 13.1607033 1.2184245 3.289774 0.3703672 NA NA
ASV42 164.2897896 -0.1230981 1.760244 -0.0699324 0.9442475 0.9919616
ASV43 96.3928250 0.3648341 3.286129 0.1110225 0.9115985 0.9919616
ASV44 154.1188333 -0.3650474 2.009368 -0.1816727 0.8558396 0.9919616
ASV45 191.8091766 -1.2744567 1.998382 -0.6377442 0.5236402 0.9919616
ASV46 74.9889042 2.1514907 3.286245 0.6546957 0.5126637 0.9919616
ASV47 52.8249013 0.9780768 3.286553 0.2975996 NA NA
ASV48 21.0604389 -21.6189123 3.390043 -6.3771802 0.0000000 0.0000000
ASV49 146.8648450 3.0450581 2.889235 1.0539323 0.2919139 0.9517866
ASV50 136.5104373 -0.3760083 2.153932 -0.1745683 0.8614189 0.9919616
ASV51 117.6774520 0.6283688 2.589345 0.2426748 0.8082573 0.9919616
ASV52 133.4977730 -0.0250548 2.019118 -0.0124088 0.9900995 0.9952295
ASV53 147.0733955 1.3014031 3.285837 0.3960644 0.6920575 0.9919616
ASV54 4.6224055 -5.1142311 3.392380 -1.5075643 0.1316661 0.8239591
ASV55 52.8793226 1.1494281 3.286533 0.3497388 NA NA
ASV56 104.9625328 -0.3409819 2.289277 -0.1489474 0.8815951 0.9919616
ASV57 108.5927540 -0.1162813 2.274199 -0.0511307 0.9592214 0.9919616
ASV58 73.4771208 -2.5146396 2.530601 -0.9936925 NA NA
ASV59 154.7492050 -1.6283106 2.134479 -0.7628610 0.4455463 0.9919616
ASV60 136.6390867 -1.1831715 2.025412 -0.5841634 0.5591104 0.9919616
ASV61 20.9664993 -21.6141366 3.390046 -6.3757658 0.0000000 0.0000000
ASV62 55.9784879 -1.6604781 3.288210 -0.5049793 NA NA
ASV63 26.8069736 25.5388291 3.323020 7.6854285 NA NA
ASV64 161.7641332 1.5093733 2.403896 0.6278862 0.5300785 0.9919616
ASV65 81.6119418 0.1705516 2.576779 0.0661879 0.9472282 0.9919616
ASV66 105.2490022 0.1401011 2.433482 0.0575723 0.9540893 0.9919616
ASV67 74.9196349 -2.8076000 2.369384 -1.1849491 NA NA
ASV68 26.8957111 5.7690474 3.025839 1.9065941 NA NA
ASV69 117.3571525 -0.1675072 2.155159 -0.0777238 0.9380477 0.9919616
ASV70 12.5135518 2.4985861 3.290537 0.7593247 NA NA
ASV71 13.4673147 4.1309178 3.294256 1.2539759 NA NA
ASV72 92.3464960 1.0178174 2.001669 0.5084845 0.6111136 0.9919616
ASV73 101.0965922 0.8143772 2.287125 0.3560702 0.7217880 0.9919616
ASV74 129.7373560 -0.3492474 1.747405 -0.1998664 0.8415851 0.9919616
ASV75 134.3798926 -0.4867884 1.749176 -0.2782959 0.7807852 0.9919616
ASV76 122.7898792 -2.0713988 1.992039 -1.0398386 0.2984149 0.9568046
ASV77 37.8789280 -0.8495222 3.288146 -0.2583591 NA NA
ASV78 138.2753466 -24.2427804 3.338224 -7.2621789 0.0000000 0.0000000
ASV79 34.5747689 2.0240758 3.287141 0.6157557 0.5380558 0.9919616
ASV80 117.7473657 1.3424219 3.285932 0.4085361 0.6828801 0.9919616
ASV81 105.9333316 0.2802472 1.885028 0.1486700 0.8818140 0.9919616
ASV82 82.2507713 1.2290321 2.395312 0.5130989 0.6078822 0.9919616
ASV83 46.6562218 -1.9613577 3.289373 -0.5962709 NA NA
ASV84 93.4199713 0.6856970 2.585498 0.2652089 0.7908485 0.9919616
ASV85 36.6760082 1.6387017 3.286995 0.4985410 NA NA
ASV86 86.0773333 0.2100168 2.419947 0.0867857 0.9308419 0.9919616
ASV87 90.7172497 -0.8728382 2.268956 -0.3846872 0.7004692 0.9919616
ASV88 32.9963492 1.9926359 3.287215 0.6061776 0.5443968 0.9919616
ASV89 102.9145267 -2.3340245 2.255945 -1.0346103 0.3008509 0.9568046
ASV90 78.3930056 0.7356288 1.958344 0.3756383 0.7071858 0.9919616
ASV91 111.0160229 0.9400974 3.285978 0.2860937 NA NA
ASV92 123.6839801 -1.9560215 2.421463 -0.8077851 0.4192143 0.9919616
ASV93 108.1357923 -1.7595728 2.262041 -0.7778695 0.4366460 0.9919616
ASV94 102.5765207 0.5920580 2.015561 0.2937435 0.7689539 0.9919616
ASV95 18.9632444 -21.4660581 3.390115 -6.3319549 0.0000000 0.0000000
ASV96 109.5682686 -0.4250654 1.749458 -0.2429698 0.8080288 0.9919616
ASV97 94.5684637 0.1747212 1.759452 0.0993043 0.9208966 0.9919616
ASV98 130.1495422 -1.3724365 2.281659 -0.6015083 0.5475015 0.9919616
ASV99 117.7058418 1.6244231 3.285932 0.4943569 0.6210542 0.9919616
ASV100 102.9937411 -2.1619775 1.962247 -1.1017869 0.2705543 0.9181599
ASV101 37.4555346 0.1752158 2.825810 0.0620055 NA NA
ASV102 30.3563940 0.8641839 3.236315 0.2670272 0.7894482 0.9919616
ASV103 104.7422520 -0.0628827 2.145861 -0.0293042 0.9766220 0.9919616
ASV104 77.2801867 1.0011760 1.953669 0.5124594 0.6083295 0.9919616
ASV105 112.9382241 -23.9538555 3.303113 -7.2519033 0.0000000 0.0000000
ASV106 116.4831763 -3.9291606 2.392126 -1.6425394 0.1004783 0.8011152
ASV107 29.7158086 -0.7888388 3.254288 -0.2423998 NA NA
ASV108 27.2525697 1.3666612 3.287526 0.4157112 NA NA
ASV109 8.5787506 2.9838815 3.293954 0.9058663 NA NA
ASV110 108.7714319 -23.9037073 3.290335 -7.2648249 0.0000000 0.0000000
ASV111 50.5524922 3.6920541 3.287360 1.1231060 NA NA
ASV112 89.3307534 -1.2525755 3.199650 -0.3914726 0.6954479 0.9919616
ASV113 61.5370469 1.2055573 2.538284 0.4749498 0.6348228 0.9919616
ASV114 107.2245667 -1.3185745 2.437291 -0.5410000 0.5885076 0.9919616
ASV115 83.7157143 -1.8099497 2.096484 -0.8633262 0.3879581 0.9919616
ASV116 15.0144655 24.7446523 3.323423 7.4455326 NA NA
ASV117 89.4470181 -0.1933697 1.740902 -0.1110744 0.9115573 0.9919616
ASV118 36.4658140 2.0453866 3.287059 0.6222544 0.5337746 0.9919616
ASV119 96.7864149 -1.3176077 2.586825 -0.5093533 0.6105046 0.9919616
ASV120 91.3182811 0.5654276 1.983503 0.2850651 0.7755942 0.9919616
ASV121 29.7605845 0.5276219 2.813538 0.1875297 0.8512454 0.9919616
ASV122 59.7716010 0.7052857 2.086231 0.3380669 0.7353127 0.9919616
ASV123 70.0608725 -0.2748272 2.397299 -0.1146403 0.9087302 0.9919616
ASV124 37.3398965 1.9666688 3.287006 0.5983162 NA NA
ASV125 72.4822695 0.0395099 3.286442 0.0120221 NA NA
ASV126 85.8037812 0.8863308 3.286138 0.2697181 0.7873771 0.9919616
ASV127 81.2073174 -3.2216932 2.207768 -1.4592536 0.1444953 0.8239591
ASV128 14.0029345 -21.0527886 3.390374 -6.2095766 0.0000000 0.0000000
ASV129 75.4741929 0.5466655 2.408316 0.2269907 0.8204310 0.9919616
ASV130 84.6584330 -1.2258763 2.825675 -0.4338349 0.6644083 0.9919616
ASV131 96.1575862 1.2092372 3.222729 0.3752215 0.7074957 0.9919616
ASV132 60.4910383 -0.1026891 2.721032 -0.0377390 0.9698958 0.9919616
ASV133 77.1898867 -0.3008719 1.761643 -0.1707904 0.8643885 0.9919616
ASV134 91.4972903 -1.7801029 2.577559 -0.6906158 0.4898070 0.9919616
ASV135 63.6047328 0.7157825 1.963856 0.3644782 0.7155010 0.9919616
ASV136 57.8012037 0.8251714 2.565021 0.3217017 0.7476787 0.9919616
ASV137 102.0358370 -2.9711336 2.241957 -1.3252414 0.1850911 0.8254637
ASV138 6.5626060 2.0304048 3.294337 0.6163318 NA NA
ASV139 76.8622643 -1.5698741 2.237371 -0.7016602 0.4828911 0.9919616
ASV140 6.2522866 -5.5479189 3.391603 -1.6357809 0.1018855 0.8011152
ASV141 81.6726360 0.5822274 2.425495 0.2400448 0.8102955 0.9919616
ASV142 55.4562266 1.1084292 2.386938 0.4643729 0.6423806 0.9919616
ASV143 90.5791912 -23.6475421 2.798983 -8.4486197 0.0000000 0.0000000
ASV144 82.3717268 -1.4105106 2.400063 -0.5876973 0.5567355 0.9919616
ASV145 35.1470250 1.7299924 3.287069 0.5263024 NA NA
ASV146 33.1813171 2.1875452 3.201299 0.6833304 0.4943981 0.9919616
ASV147 76.9935825 -9.1714264 2.113207 -4.3400514 0.0000142 0.0003251
ASV148 47.4419723 0.4534362 2.245278 0.2019510 0.8399550 0.9919616
ASV149 26.7691252 -0.3444786 2.826824 -0.1218606 0.9030094 0.9919616
ASV150 68.3408444 0.5461074 2.555630 0.2136880 0.8307904 0.9919616
ASV151 69.3597547 0.6442583 2.733715 0.2356713 0.8136877 0.9919616
ASV152 22.9817179 0.7227008 3.288079 0.2197942 NA NA
ASV153 67.3421817 -1.4197520 2.735684 -0.5189751 0.6037781 0.9919616
ASV154 33.1585356 3.3676116 3.287988 1.0242165 NA NA
ASV155 58.9457682 0.8087817 2.405926 0.3361623 0.7367485 0.9919616
ASV156 80.0176808 2.9913610 2.494061 1.1993936 0.2303749 0.8845429
ASV157 68.6227101 0.0878385 2.404225 0.0365351 0.9708557 0.9919616
ASV158 62.2145850 0.6289519 2.238705 0.2809445 0.7787529 0.9919616
ASV159 77.5779727 1.0698959 3.173304 0.3371552 0.7359999 0.9919616
ASV160 48.8186867 -1.4597874 2.391080 -0.6105140 0.5415214 0.9919616
ASV161 23.9642660 0.6512551 3.288006 0.1980699 NA NA
ASV162 17.1081541 1.0647190 3.288816 0.3237394 NA NA
ASV163 12.4978267 2.3706809 3.290406 0.7204829 0.4712278 0.9919616
ASV164 46.7524583 0.2694123 3.286885 0.0819659 0.9346739 0.9919616
ASV165 26.2865534 -20.7220813 3.389912 -6.1128676 0.0000000 0.0000000
ASV166 57.2151267 0.5838973 2.096336 0.2785324 0.7806037 0.9919616
ASV167 82.8643127 -2.9766661 2.534423 -1.1744947 0.2401969 0.8845429
ASV168 79.9113593 -1.1468395 2.282712 -0.5024022 0.6153846 0.9919616
ASV169 30.5030962 2.4961217 2.952891 0.8453145 0.3979353 0.9919616
ASV170 28.6949520 2.4676398 3.287656 0.7505773 0.4529071 0.9919616
ASV171 65.7534356 -8.9438292 2.569298 -3.4810399 0.0004995 0.0107664
ASV172 65.3240072 3.8085358 3.287009 1.1586630 0.2465936 0.8941898
ASV173 27.8082679 -21.4318902 3.389883 -6.3223099 0.0000000 0.0000000
ASV174 57.5447541 0.6523279 2.555788 0.2552355 0.7985412 0.9919616
ASV175 57.5698035 -1.9116103 2.226085 -0.8587318 0.3904885 0.9919616
ASV176 31.2072349 2.5723920 3.207363 0.8020271 NA NA
ASV177 58.1177857 0.0114375 2.261235 0.0050581 0.9959643 0.9968603
ASV178 61.4946545 0.2436578 2.571961 0.0947362 0.9245244 0.9919616
ASV179 50.0383825 26.3335235 3.322781 7.9251453 NA NA
ASV180 19.9949146 -0.6167675 3.025000 -0.2038901 NA NA
ASV181 21.3556423 1.6635071 3.288104 0.5059168 0.6129151 0.9919616
ASV182 23.8275322 1.9593435 3.241651 0.6044277 0.5455594 0.9919616
ASV183 40.2638844 -1.9626797 2.740299 -0.7162284 0.4738504 0.9919616
ASV184 53.7619259 1.0205699 3.286528 0.3105313 0.7561570 0.9919616
ASV185 22.3103569 3.5452616 3.289507 1.0777487 NA NA
ASV186 53.4372593 -0.4119162 2.732149 -0.1507664 0.8801600 0.9919616
ASV187 18.1301853 2.8595457 3.289325 0.8693411 NA NA
ASV188 7.9598140 -5.8981147 3.391125 -1.7392797 0.0819856 0.8011152
ASV189 53.8165630 1.1481840 2.550386 0.4502001 0.6525662 0.9919616
ASV190 70.8269217 -1.4271641 2.398874 -0.5949309 0.5518896 0.9919616
ASV191 56.6177617 -0.5870391 2.566645 -0.2287185 0.8190877 0.9919616
ASV192 46.6467565 1.0898651 2.397980 0.4544930 0.6494740 0.9919616
ASV193 48.7987065 4.9216474 3.289188 1.4963108 NA NA
ASV194 40.7746917 1.6752866 2.573610 0.6509481 0.5150800 0.9919616
ASV195 48.1589724 1.2649562 3.055738 0.4139610 0.6789027 0.9919616
ASV196 4.3287634 6.4862999 3.325683 1.9503662 NA NA
ASV197 54.5010051 1.2348629 3.286494 0.3757386 0.7071112 0.9919616
ASV198 50.8565292 -22.8484206 3.389657 -6.7406292 0.0000000 0.0000000
ASV199 3.8264294 6.2027020 3.326226 1.8647865 NA NA
ASV200 47.7774942 0.3886691 3.072326 0.1265064 0.8993311 0.9919616
ASV201 44.2030046 -1.5236868 2.560120 -0.5951623 0.5517350 0.9919616
ASV202 22.6950177 0.4795746 3.288250 0.1458449 NA NA
ASV203 54.1301506 -1.3265726 2.723668 -0.4870538 0.6262202 0.9919616
ASV204 48.8049748 1.8067306 3.286622 0.5497226 0.5825096 0.9919616
ASV205 55.2765466 -3.5158315 2.488553 -1.4128016 0.1577141 0.8239591
ASV206 26.3535301 2.2435286 3.287745 0.6823914 0.4949915 0.9919616
ASV207 38.4959885 1.0155406 3.286956 0.3089608 0.7573513 0.9919616
ASV208 14.9229529 0.4113114 3.255043 0.1263613 NA NA
ASV209 11.5976292 2.0043164 3.290469 0.6091279 NA NA
ASV210 18.3863950 1.2955764 3.254848 0.3980452 NA NA
ASV211 29.9647993 0.2875438 3.025092 0.0950529 0.9242728 0.9919616
ASV212 34.0261320 -1.6277164 3.230154 -0.5039129 0.6143226 0.9919616
ASV213 46.1236192 3.1199158 3.287113 0.9491356 0.3425517 0.9919616
ASV214 15.6301085 3.3025154 3.290695 1.0035921 NA NA
ASV215 46.9193516 -8.4568881 2.448534 -3.4538574 0.0005526 0.0112853
ASV216 34.4651476 0.3815038 3.210697 0.1188227 0.9054158 0.9919616
ASV217 53.9456636 -3.4929526 2.198408 -1.5888553 0.1120931 0.8119757
ASV218 5.5393251 -5.3755306 3.391884 -1.5848214 0.1130069 0.8119757
ASV219 35.0962495 -0.0550694 3.287560 -0.0167508 0.9866354 0.9952295
ASV220 28.2677789 0.4431015 3.050036 0.1452775 NA NA
ASV221 52.1819645 1.1771428 3.286545 0.3581703 0.7202159 0.9919616
ASV222 22.5137246 2.2530047 3.216704 0.7004078 0.4836727 0.9919616
ASV223 48.3942311 -0.7778237 2.123283 -0.3663306 0.7141184 0.9919616
ASV224 23.3224540 -0.5419084 3.273306 -0.1655538 0.8685081 0.9919616
ASV225 17.1820961 3.2306344 3.290092 0.9819283 0.3261352 0.9919616
ASV226 13.9613476 1.4794383 3.289487 0.4497474 NA NA
ASV227 15.4459582 0.7451223 3.253907 0.2289931 NA NA
ASV228 18.0950799 -21.4099196 3.390150 -6.3153305 0.0000000 0.0000000
ASV229 8.9150563 1.2036123 3.291834 0.3656358 NA NA
ASV230 18.1872632 2.2421373 3.288772 0.6817553 0.4953937 0.9919616
ASV231 36.9628833 1.2164316 3.286991 0.3700745 0.7113270 0.9919616
ASV232 8.4767189 1.2217978 3.292156 0.3711239 NA NA
ASV233 6.3483555 -5.5716495 3.391567 -1.6427952 0.1004253 0.8011152
ASV234 14.2772380 1.1927770 3.289440 0.3626080 NA NA
ASV235 24.2052086 0.7532608 3.013809 0.2499365 0.8026365 0.9919616
ASV236 9.2810667 24.0945809 3.323988 7.2486960 NA NA
ASV237 36.6586500 -0.1201238 2.605021 -0.0461124 0.9632206 0.9919616
ASV238 16.3895841 2.5984305 3.289433 0.7899326 0.4295671 0.9919616
ASV239 27.5059995 1.9752325 3.207687 0.6157809 0.5380391 0.9919616
ASV240 41.7755558 -0.3194602 2.588633 -0.1234088 0.9017834 0.9919616
ASV241 37.3896336 0.5503419 2.790097 0.1972483 0.8436332 0.9919616
ASV242 9.4064437 1.0695397 3.165097 0.3379169 0.7354258 0.9919616
ASV243 21.2411784 2.5452159 3.288484 0.7739784 0.4389435 0.9919616
ASV244 49.5806522 0.0121312 3.082893 0.0039350 0.9968603 0.9968603
ASV245 10.6788339 0.6902524 3.291137 0.2097306 NA NA
ASV246 13.0591888 0.0089209 3.290996 0.0027107 NA NA
ASV247 9.9327054 -6.2175365 3.390779 -1.8336600 0.0667045 0.8011152
ASV248 32.5375699 3.7352951 3.288473 1.1358753 NA NA
ASV249 35.4814086 -0.8575129 2.796340 -0.3066554 0.7591057 0.9919616
ASV250 52.1309124 0.8530368 3.286586 0.2595510 0.7952101 0.9919616
ASV251 19.5520259 2.0659194 3.288457 0.6282338 0.5298508 0.9919616
ASV252 39.2518547 -0.5221785 2.986130 -0.1748679 0.8611834 0.9919616
ASV253 43.6216423 2.6592966 3.286972 0.8090414 0.4184913 0.9919616
ASV254 25.4158867 1.2303131 3.021052 0.4072465 0.6838269 0.9919616
ASV255 12.4987439 24.4982886 3.323607 7.3709940 NA NA
ASV256 24.2134821 -1.5551125 3.010744 -0.5165209 0.6054906 0.9919616
ASV257 39.0395996 -1.9672033 3.290153 -0.5979064 0.5499024 0.9919616
ASV258 27.0924502 0.6565059 2.833005 0.2317348 0.8167440 0.9919616
ASV259 35.6569182 -8.0610652 2.652187 -3.0394032 0.0023705 0.0459872
ASV260 16.9292489 0.0401853 3.289687 0.0122155 NA NA
ASV261 19.3507345 1.7627943 3.288394 0.5360654 NA NA
ASV262 34.2345214 1.4738692 3.190006 0.4620270 0.6440619 0.9919616
ASV263 29.8613336 0.6897162 3.287487 0.2098004 0.8338234 0.9919616
ASV264 32.6662047 -2.1561104 3.201494 -0.6734701 0.5006482 0.9919616
ASV265 7.6861822 23.8374323 3.324296 7.1706707 NA NA
ASV266 14.8354575 0.7898259 3.289476 0.2401069 0.8102474 0.9919616
ASV267 31.9668768 -0.5587962 3.288257 -0.1699369 0.8650598 0.9919616
ASV268 31.9146218 -0.8801095 3.238485 -0.2717658 0.7858021 0.9919616
ASV269 11.4759410 4.5683018 3.298466 1.3849778 0.1660593 0.8239591
ASV270 20.4304992 1.6167745 3.288219 0.4916870 0.6229407 0.9919616
ASV271 4.6138978 -5.1113368 3.392386 -1.5067084 0.1318854 0.8239591
ASV272 7.5106329 1.2960818 3.292987 0.3935885 NA NA
ASV273 20.1399236 1.9597788 3.288329 0.5959802 0.5511885 0.9919616
ASV274 25.6686758 -0.2381852 3.260728 -0.0730466 0.9417690 0.9919616
ASV275 21.9567158 1.6583223 3.288030 0.5043513 0.6140145 0.9919616
ASV276 21.5756826 -21.6536144 3.390027 -6.3874459 0.0000000 0.0000000
ASV277 34.6437893 2.3887362 2.953981 0.8086498 0.4187166 0.9919616
ASV278 45.0953610 -1.0949873 2.758066 -0.3970127 0.6913581 0.9919616
ASV279 39.5523689 -0.7212364 3.287885 -0.2193618 0.8263682 0.9919616
ASV280 31.4389091 0.7801898 3.235443 0.2411385 NA NA
ASV281 16.6110098 1.7512755 3.288877 0.5324843 0.5943906 0.9919616
ASV282 30.2849348 0.9869484 3.287371 0.3002243 0.7640061 0.9919616
ASV283 7.4610846 1.0107439 3.293197 0.3069187 0.7589052 0.9919616
ASV284 22.4614687 0.6796489 3.256092 0.2087315 0.8346578 0.9919616
ASV285 10.2558500 -0.4690832 3.293867 -0.1424111 NA NA
ASV286 33.6273193 -0.0593059 3.235164 -0.0183316 0.9853743 0.9952295
ASV287 10.9558645 3.4559180 3.293399 1.0493469 NA NA
ASV288 23.1469912 -0.9516819 3.252442 -0.2926054 0.7698238 0.9919616
ASV289 33.7411213 0.5850212 3.287291 0.1779645 0.8587509 0.9919616
ASV290 5.0439992 -5.2372063 3.392136 -1.5439260 0.1226063 0.8239591
ASV291 17.7192392 0.2909657 3.289209 0.0884607 0.9295105 0.9919616
ASV292 16.3203729 1.6685669 3.288923 0.5073293 0.6119238 0.9919616
ASV293 21.7574716 -0.5882902 3.263430 -0.1802674 0.8569426 0.9919616
ASV294 16.4221681 2.7611878 3.289603 0.8393681 0.4012628 0.9919616
ASV295 6.0180935 -5.4934464 3.391688 -1.6196791 0.1053012 0.8011152
ASV296 7.8064369 3.9845315 3.299557 1.2075958 NA NA
ASV297 9.9510813 24.1891958 3.323889 7.2773777 NA NA
ASV298 1.4781120 4.9375253 3.331789 1.4819440 NA NA
ASV299 14.9884137 -21.1466952 3.369857 -6.2752507 0.0000000 0.0000000
ASV300 31.0897454 1.4297142 3.287267 0.4349249 0.6636170 0.9919616
ASV301 7.9523599 23.8837630 3.324237 7.1847365 NA NA
ASV302 18.2433368 -0.8310908 3.290990 -0.2525352 0.8006274 0.9919616
ASV303 2.0995546 -3.9849366 3.395931 -1.1734444 0.2406176 0.8845429
ASV304 23.9470950 1.2435789 3.287824 0.3782377 NA NA
ASV305 19.4090722 2.0256975 3.288461 0.6160018 0.5378934 0.9919616
ASV306 6.8574877 1.2261718 3.293736 0.3722739 NA NA
ASV307 25.5499679 1.5054092 3.025120 0.4976362 0.6187405 0.9919616
ASV308 15.3383733 0.1406839 2.996615 0.0469476 0.9625550 0.9919616
ASV309 20.5325467 -4.2107017 2.789221 -1.5096336 0.1311369 0.8239591
ASV310 12.6632098 -0.2869019 2.970557 -0.0965818 NA NA
ASV311 6.2437627 -5.5469686 3.391604 -1.6354999 0.1019443 0.8011152
ASV312 10.0826985 1.8543329 3.291134 0.5634328 NA NA
ASV313 22.3408421 -0.0531274 3.288762 -0.0161542 NA NA
ASV314 7.8861252 2.6442201 3.293810 0.8027847 NA NA
ASV315 15.1683373 1.9318736 3.289260 0.5873277 NA NA
ASV316 13.3758267 3.6666508 3.292579 1.1136104 0.2654464 0.9179515
ASV317 8.4905960 1.2789185 3.292129 0.3884777 NA NA
ASV318 30.7943441 0.6231298 3.239203 0.1923714 0.8474513 0.9919616
ASV319 28.2112952 -0.2694610 2.848642 -0.0945928 NA NA
ASV320 18.4037310 2.0780666 3.288649 0.6318907 NA NA
ASV321 22.5055690 2.6995797 3.288432 0.8209323 0.4116848 0.9919616
ASV322 14.8047208 -2.2036178 3.244258 -0.6792363 0.4969881 0.9919616
ASV323 9.4589579 0.8903929 3.291664 0.2704993 NA NA
ASV324 21.0513109 1.6556005 3.229021 0.5127252 0.6081435 0.9919616
ASV325 6.4577070 1.6873393 3.294222 0.5122118 NA NA
ASV326 22.5854420 0.0478903 3.252817 0.0147227 0.9882534 0.9952295
ASV327 19.9823089 -2.8710827 2.984753 -0.9619164 0.3360916 0.9919616
ASV328 19.2104796 -1.3709899 2.624961 -0.5222896 0.6014687 0.9919616
ASV329 5.0530925 6.7096768 3.325228 2.0178098 0.0436111 0.7356998
ASV330 25.7234206 -21.8979274 3.388500 -6.4624244 0.0000000 0.0000000
ASV331 7.2275408 4.6623864 3.307107 1.4098081 0.1585964 0.8239591
ASV332 15.4107524 1.7090247 3.289137 0.5195966 0.6033448 0.9919616
ASV333 9.8369510 2.2011020 3.291550 0.6687129 NA NA
ASV334 31.8099597 -22.1947664 3.351546 -6.6222471 0.0000000 0.0000000
ASV335 27.1051847 2.0826697 3.229479 0.6448934 0.5189963 0.9919616
ASV336 6.5372984 -5.6122495 3.391506 -1.6547956 0.0979660 0.8011152
ASV337 2.5298948 -4.2431824 3.394860 -1.2498844 0.2113418 0.8814503
ASV338 9.6617263 3.5721724 3.294911 1.0841482 0.2782991 0.9229064
ASV339 25.5256802 0.0969723 3.288212 0.0294909 0.9764731 0.9919616
ASV340 3.2319751 6.0645677 3.326761 1.8229647 0.0683087 0.8011152
ASV341 22.3250017 1.1194853 3.256602 0.3437587 0.7310277 0.9919616
ASV342 6.8223288 0.7725510 3.294216 0.2345174 NA NA
ASV343 2.2755830 5.5601056 3.328540 1.6704339 NA NA
ASV344 11.7327809 0.5186811 3.290811 0.1576150 0.8747602 0.9919616
ASV345 10.3763412 24.2051116 3.323832 7.2822900 NA NA
ASV346 13.5946642 0.5370045 3.290062 0.1632202 0.8703451 0.9919616
ASV347 3.5719895 22.7924678 3.326353 6.8520890 NA NA
ASV348 10.9576087 1.1924889 3.290646 0.3623875 0.7170625 0.9919616
ASV349 12.6577918 2.3155956 3.290286 0.7037673 0.4815777 0.9919616
ASV350 25.2809479 0.5876692 3.261714 0.1801719 0.8570176 0.9919616
ASV351 1.0759840 -3.0155157 3.402181 -0.8863478 0.3754301 0.9919616
ASV352 2.5994486 -4.2911586 3.394681 -1.2640830 0.2062002 0.8791834
ASV353 8.4659340 -0.2966824 3.294978 -0.0900408 NA NA
ASV354 6.6175349 -5.6308267 3.391479 -1.6602865 0.0968568 0.8011152
ASV355 9.6891171 0.5318415 3.291923 0.1615595 NA NA
ASV356 16.1466007 2.4615437 3.289362 0.7483348 NA NA
ASV357 1.1744862 1.4962560 3.333397 0.4488682 0.6535267 0.9919616
ASV358 4.5719080 -5.0969026 3.392416 -1.5024402 0.1329835 0.8239591
ASV359 4.1099175 -4.9429479 3.392757 -1.4569117 0.1451408 0.8239591
ASV360 9.1161164 0.2616687 3.292802 0.0794669 NA NA
ASV361 2.9150425 5.9160224 3.327222 1.7780667 NA NA
ASV362 10.5798567 2.2664968 3.291188 0.6886562 0.4910397 0.9919616
ASV363 8.1041759 3.5559606 2.750697 1.2927492 0.1960978 0.8646131
ASV364 4.3866723 -5.0390297 3.392540 -1.4853264 0.1374574 0.8239591
ASV365 17.5116212 2.8908571 3.289501 0.8788132 NA NA
ASV366 11.1577084 -0.1605420 3.292323 -0.0487625 NA NA
ASV367 15.7930769 0.1635221 2.781788 0.0587831 0.9531249 0.9919616
ASV368 22.6282192 2.7141937 3.288428 0.8253773 0.4091574 0.9919616
ASV369 5.7247150 1.3677967 3.295281 0.4150774 NA NA
ASV370 1.1381014 -3.0806698 3.401618 -0.9056484 0.3651220 0.9919616
ASV371 5.1695403 1.9520955 3.296663 0.5921428 NA NA
ASV372 4.4594402 23.0867014 3.325592 6.9421330 NA NA
ASV373 2.1918493 -4.0429070 3.395674 -1.1906053 0.2338086 0.8845429
ASV374 15.9080229 1.6743958 3.289014 0.5090875 0.6106909 0.9919616
ASV375 2.3550960 -4.1492944 3.395227 -1.2220962 0.2216713 0.8845429
ASV376 6.9135138 -5.6946465 3.391388 -1.6791490 0.0931230 0.8011152
ASV377 8.1874175 1.4892460 3.292338 0.4523369 NA NA
ASV378 2.4720810 5.6770332 3.328070 1.7058032 NA NA
ASV379 9.7603254 3.7740186 3.295689 1.1451381 NA NA
ASV380 3.7848130 6.2923456 3.326140 1.8917862 NA NA
ASV381 8.7680420 2.1490889 3.292217 0.6527786 NA NA
ASV382 9.7940068 3.4540112 3.294331 1.0484712 NA NA
ASV383 23.3290658 -7.4488281 3.143582 -2.3695349 0.0178105 0.3290697
ASV384 8.8949448 1.6494805 3.291812 0.5010859 NA NA
ASV385 3.1147029 6.0115655 3.326920 1.8069463 NA NA
ASV386 6.6812132 -5.6445462 3.391459 -1.6643415 0.0960442 0.8011152
ASV387 6.4190207 0.1818046 3.296135 0.0551569 NA NA
ASV388 6.4079418 -5.5823730 3.391551 -1.6459649 0.0997710 0.8011152
ASV389 11.3484217 2.7554524 3.291442 0.8371566 0.4025046 0.9919616
ASV390 7.0834173 1.6255554 3.293411 0.4935780 0.6216042 0.9919616
ASV391 6.8358996 -5.6774673 3.391412 -1.6740716 0.0941165 0.8011152
ASV392 3.1347386 -4.5518680 3.393806 -1.3412279 0.1798465 0.8239591
ASV393 10.2962006 -6.2675505 3.390732 -1.8484359 0.0645393 0.8011152
ASV394 3.4740917 6.1690542 3.326464 1.8545382 NA NA
ASV395 10.2387843 0.5422407 3.291565 0.1647364 NA NA
ASV396 23.2932545 1.2940604 3.287885 0.3935845 0.6938879 0.9919616
ASV397 4.3745643 2.4271239 3.299739 0.7355502 NA NA
ASV398 4.7484277 5.4350098 3.324544 1.6348136 NA NA
ASV399 3.5911397 -4.7489913 3.393242 -1.3995438 0.1616500 0.8239591
ASV400 13.8403938 1.3775190 3.289532 0.4187583 0.6753928 0.9919616
ASV401 7.1472598 1.7943223 3.293413 0.5448216 NA NA
ASV402 6.0998737 -5.5128067 3.391657 -1.6254020 0.1040769 0.8011152
ASV403 7.1915050 -1.5644982 3.167674 -0.4938950 NA NA
ASV404 14.0543729 0.9313149 3.289614 0.2831077 0.7770943 0.9919616
ASV405 13.8193594 -6.6931086 3.358081 -1.9931350 0.0462467 0.7476544
ASV406 3.3270091 -1.4441478 3.326308 -0.4341594 0.6641727 0.9919616
ASV407 9.7024443 1.6426026 3.291280 0.4990771 0.6177251 0.9919616
ASV408 11.7020946 -0.6230667 3.293328 -0.1891906 0.8499434 0.9919616
ASV409 8.9737252 1.9195733 3.291875 0.5831247 0.5598094 0.9919616
ASV410 2.7476842 5.8320359 3.327504 1.7526758 NA NA
ASV411 9.7297698 0.1064829 3.292665 0.0323394 0.9742014 0.9919616
ASV412 6.9994500 -5.7105018 3.391367 -1.6838351 0.0922135 0.8011152
ASV413 0.0000000 0.0000000 0.000000 0.0000000 1.0000000 NA
ASV414 4.4228487 1.0500917 3.298443 0.3183598 NA NA
ASV415 2.5690254 5.7325532 3.327861 1.7225941 NA NA
ASV416 12.3453881 -0.6414907 3.292975 -0.1948058 0.8455449 0.9919616
ASV417 6.6150192 -5.6323120 3.391477 -1.6607255 0.0967686 0.8011152
ASV418 3.8036444 -4.8298986 3.393032 -1.4234757 0.1545983 0.8239591
ASV419 16.0023439 2.6099305 3.289538 0.7934033 0.4275429 0.9919616
ASV420 2.4239814 5.6495918 3.328177 1.6975033 NA NA
ASV421 2.3881343 -4.1661621 3.395160 -1.2270887 0.2197892 0.8845429
ASV422 6.7740307 -5.6656290 3.391429 -1.6705727 0.0948061 0.8011152
ASV423 7.6803861 -5.8441031 3.391191 -1.7233186 0.0848309 0.8011152
ASV424 4.3557170 0.9889644 3.298775 0.2997975 0.7643316 0.9919616
ASV425 5.7186088 6.8864402 3.324914 2.0711635 0.0383435 0.6762403
ASV426 8.5202847 3.6565795 3.296575 1.1092055 0.2673415 0.9179515
ASV427 2.5642318 0.3011807 3.311220 0.0909576 NA NA
ASV428 8.2709097 -5.9519622 3.391061 -1.7551916 0.0792265 0.8011152
ASV429 5.3842392 2.7507337 3.298058 0.8340464 NA NA
ASV430 3.3595531 6.1204137 3.326599 1.8398409 NA NA
ASV431 13.5419296 2.0505148 3.289780 0.6232985 0.5330884 0.9919616
ASV432 2.1745387 -4.0352213 3.395707 -1.1883301 0.2347034 0.8845429
ASV433 3.4265072 -4.6755238 3.393444 -1.3778110 0.1682617 0.8239591
ASV434 2.1011621 -3.9762197 3.395971 -1.1708639 0.2416535 0.8845429
ASV435 1.9410121 -3.8718021 3.396464 -1.1399508 0.2543068 0.9136207
ASV436 6.1651223 4.6268684 3.310374 1.3976874 NA NA
ASV437 7.9489448 0.7476275 3.293017 0.2270342 0.8203971 0.9919616
ASV438 3.9901374 2.2086291 3.300447 0.6691909 NA NA
ASV439 1.4523470 -3.4419637 3.398915 -1.0126654 0.3112200 0.9759373
ASV440 1.1979627 4.6338383 3.333960 1.3898904 NA NA
ASV441 3.0751319 -4.5297205 3.393875 -1.3346752 0.1819827 0.8239591
ASV442 7.7704493 -0.5019681 3.296708 -0.1522634 0.8789792 0.9919616
ASV443 3.1325607 6.0167761 3.326904 1.8085211 NA NA
ASV444 7.8490268 1.7870771 3.292708 0.5427377 NA NA
ASV445 0.4196783 -1.7327632 3.419279 -0.5067628 0.6123213 0.9919616
ASV446 11.1654858 -6.3848205 3.390627 -1.8830798 0.0596896 0.8011152
ASV447 6.7424884 4.1324456 3.303034 1.2511060 NA NA
ASV448 3.3273285 -4.6401090 3.393544 -1.3673341 0.1715206 0.8239591
ASV449 3.7199159 6.2651016 3.326209 1.8835562 NA NA
ASV450 4.5030560 6.5411053 3.325565 1.9669157 NA NA
ASV451 2.2228125 -4.0566065 3.395614 -1.1946605 0.2322197 0.8845429
ASV452 7.1174911 1.0877348 3.293514 0.3302657 NA NA
ASV453 8.1716461 2.9002726 3.294150 0.8804314 0.3786257 0.9919616
ASV454 10.3049731 0.8196009 3.291217 0.2490267 NA NA
ASV455 7.3939638 3.2358032 3.296264 0.9816578 NA NA
ASV456 0.7229418 3.9046857 3.341463 1.1685556 NA NA
ASV457 5.3749674 2.4530743 3.297154 0.7439975 NA NA
ASV458 3.9974429 22.8884702 3.325947 6.8817911 NA NA
ASV459 2.6255649 5.7665333 3.327736 1.7328698 NA NA
ASV460 3.8859625 -4.8643999 3.392946 -1.4336803 0.1516635 0.8239591
ASV461 4.3263698 -5.0178655 3.392587 -1.4790676 0.1391222 0.8239591
ASV462 2.2113514 5.5171646 3.328722 1.6574424 NA NA
ASV463 5.6843428 -5.4116591 3.391822 -1.5955019 0.1106000 0.8119757
ASV464 15.0021654 1.3422671 3.289220 0.4080806 0.6832145 0.9919616
ASV465 7.9773806 1.2336232 3.292576 0.3746681 NA NA
ASV466 1.8739961 -3.8099752 3.396774 -1.1216453 0.2620133 0.9179515
ASV467 11.8041325 2.4503960 3.290782 0.7446244 0.4564988 0.9919616
ASV468 4.0685959 -4.9270168 3.392795 -1.4522000 0.1464460 0.8239591
ASV469 8.8361550 -6.0470282 3.390954 -1.7832821 0.0745404 0.8011152
ASV470 3.2745011 6.0834242 3.326706 1.8286632 NA NA
ASV471 12.8718050 1.8700006 3.289916 0.5684038 0.5697608 0.9919616
ASV472 5.1334321 -5.2625545 3.392088 -1.5514205 0.1208009 0.8239591
ASV473 4.7372088 3.9529295 3.308363 1.1948294 NA NA
ASV474 1.3573802 1.4857576 3.326681 0.4466186 NA NA
ASV475 0.2235279 -0.9057717 3.424764 -0.2644771 0.7914123 0.9919616
ASV476 3.0699998 -4.5230462 3.393896 -1.3327005 0.1826301 0.8239591
ASV477 2.4847749 -4.2228624 3.394937 -1.2438705 0.2135472 0.8814503
ASV478 10.2826705 1.3143547 3.290955 0.3993840 NA NA
ASV479 6.6547053 1.9096359 3.294104 0.5797134 NA NA
ASV480 10.4738435 1.7826559 3.290896 0.5416932 NA NA
ASV481 0.2249523 -1.1882658 3.423912 -0.3470492 0.7285544 0.9919616
ASV482 2.5087220 -4.2308116 3.394907 -1.2462232 0.2126825 0.8814503
ASV483 3.6658572 -4.7790901 3.393163 -1.4084471 0.1589987 0.8239591
ASV484 11.2542792 2.5356928 3.291159 0.7704558 0.4410296 0.9919616
ASV485 6.1345981 23.5259778 3.324751 7.0760117 NA NA
ASV486 9.2624884 -6.1148235 3.390882 -1.8033134 0.0713390 0.8011152
ASV487 0.8368765 -2.6543018 3.405805 -0.7793463 0.4357757 0.9919616
ASV488 4.9177500 -5.1987673 3.392210 -1.5325607 0.1253841 0.8239591
ASV489 8.1782436 0.2400018 3.293691 0.0728671 0.9419118 0.9919616
ASV490 1.7402000 5.1741333 3.330387 1.5536130 NA NA
ASV491 5.8695864 1.1652273 3.295163 0.3536175 NA NA
ASV492 4.7635678 4.0566959 3.309351 1.2258284 NA NA
ASV493 12.1468607 0.9961102 3.290228 0.3027481 0.7620818 0.9919616
ASV494 3.3296712 6.1075470 3.326636 1.8359529 0.0663646 0.8011152
ASV495 5.5991049 3.3507677 3.300292 1.0152943 NA NA
ASV496 8.3510393 2.2072796 3.292637 0.6703683 0.5026230 0.9919616
ASV497 11.3134399 5.4331437 3.307194 1.6428258 0.1004190 0.8011152
ASV498 3.3568899 -4.6533701 3.393506 -1.3712572 0.1702948 0.8239591
ASV499 7.8805435 2.4221323 3.293401 0.7354501 NA NA
ASV500 2.0580684 3.6708641 3.327078 1.1033296 NA NA
ASV501 0.1108584 1.3643242 3.418344 0.3991184 0.6898059 0.9919616
ASV502 3.3367237 2.3049513 3.303658 0.6976967 NA NA
ASV503 2.3020737 -4.1059955 3.395405 -1.2092800 0.2265553 0.8845429
ASV504 2.1998060 2.3223581 3.313194 0.7009424 NA NA
ASV505 1.5662804 5.0142549 3.331309 1.5051904 NA NA
ASV506 4.6581141 0.8176914 3.298194 0.2479210 NA NA
ASV507 2.8067152 5.8610626 3.327405 1.7614517 NA NA
ASV508 1.8814542 5.2848882 3.329806 1.5871460 NA NA
ASV509 4.7639619 -5.1581567 3.392290 -1.5205528 0.1283721 0.8239591
ASV510 0.3779434 -1.5609347 3.421622 -0.4561973 0.6482481 0.9919616
ASV511 0.6786105 3.8101904 3.342741 1.1398400 NA NA
ASV512 3.4953751 -4.7085942 3.393352 -1.3875939 0.1652607 0.8239591
ASV513 2.1844963 1.0296233 3.311992 0.3108773 NA NA
ASV514 2.2897368 5.5693694 3.328501 1.6732365 NA NA
ASV515 3.2428620 -4.6014779 3.393657 -1.3559053 0.1751293 0.8239591
ASV516 0.7765836 0.4665551 3.365866 0.1386137 0.8897554 0.9919616
ASV517 3.1752097 0.0522770 3.307868 0.0158038 NA NA
ASV518 9.8187479 3.9116049 3.296317 1.1866591 0.2353621 0.8845429
ASV519 1.3976231 4.8560517 3.332328 1.4572552 NA NA
ASV520 2.4468189 1.8987541 3.308754 0.5738578 NA NA
ASV521 1.4049190 1.5595255 3.324962 0.4690357 NA NA
ASV522 6.9852153 3.3538173 3.217753 1.0422854 0.2972794 0.9568046
ASV523 0.8317690 -2.6343256 3.406033 -0.7734292 0.4392684 0.9919616
ASV524 2.4756094 1.0398993 3.308792 0.3142837 NA NA
ASV525 2.3576093 -4.1442696 3.395248 -1.2206089 0.2222341 0.8845429
ASV526 2.6862027 -4.3273613 3.394550 -1.2747968 0.2023812 0.8791834
ASV527 5.4731279 6.4792147 3.325094 1.9485809 NA NA
ASV528 1.9096601 -3.8351181 3.396646 -1.1290897 0.2588600 0.9179515
ASV529 4.7518282 1.8717140 3.297526 0.5676116 NA NA
ASV530 2.6115207 2.4179023 3.309369 0.7306234 NA NA
ASV531 0.7152892 -2.4219293 3.408660 -0.7105224 0.4773802 0.9919616
ASV532 3.8696715 3.2796212 3.306362 0.9919124 NA NA
ASV533 0.9355717 4.2764928 3.337169 1.2814732 NA NA
ASV534 3.6924976 4.2743293 3.319008 1.2878335 NA NA
ASV535 3.1687303 1.3049795 3.303241 0.3950603 NA NA
ASV536 0.6086405 -2.2016835 3.411820 -0.6453106 0.5187259 0.9919616
ASV537 5.8014094 1.2609867 3.295224 0.3826710 0.7019637 0.9919616
ASV538 1.9427797 -0.4474780 3.329109 -0.1344137 0.8930754 0.9919616
ASV539 3.9430191 -4.8848708 3.392896 -1.4397350 0.1499424 0.8239591
ASV540 1.2851406 1.8692532 3.329799 0.5613712 NA NA
ASV541 4.4770362 5.5667952 3.325004 1.6742220 0.0940870 0.8011152
ASV542 0.6891409 -2.3415078 3.409759 -0.6867077 0.4922669 0.9919616
ASV543 0.8174920 -2.6208135 3.406189 -0.7694269 0.4416399 0.9919616
ASV544 3.2208127 -4.5886563 3.393695 -1.3521121 0.1763395 0.8239591
ASV545 1.1798472 -3.1420619 3.401110 -0.9238343 0.3555726 0.9919616
ASV546 1.2847054 1.3156473 3.329223 0.3951815 NA NA
ASV547 2.2870624 3.1179216 3.318576 0.9395361 0.3474556 0.9919616
ASV548 1.8426105 -3.7853417 3.396901 -1.1143517 0.2651284 0.9179515
ASV549 4.1550795 2.6360726 3.301294 0.7984969 0.4245822 0.9919616
ASV550 1.1482260 2.2035884 3.335215 0.6607036 NA NA
ASV551 2.6363579 -4.3017006 3.394642 -1.2672028 0.2050828 0.8791834
ASV552 0.9606208 0.9642660 3.346299 0.2881590 0.7732250 0.9919616
ASV553 0.4904125 -1.8835702 3.416901 -0.5512510 0.5814617 0.9919616
ASV554 1.6106202 0.5682425 3.323866 0.1709583 NA NA
ASV555 0.4349425 1.8073549 3.379080 0.5348661 0.5927425 0.9919616
ASV556 0.6202405 -0.1226365 3.397837 -0.0360925 0.9712086 0.9919616
ASV557 0.8806597 2.1304638 3.343906 0.6371183 0.5240478 0.9919616
#umbral
alpha = 0.05
# Ordenamos la tabla de resultados
res.PR.sCR  = res.PR.sCR [order(res.PR.sCR $padj, na.last=NA), ]
# Filtramos según nuestro umbral alpha
sigtab = res.PR.sCR [(res.PR.sCR $padj < alpha), ]
# Le agregamos la taxonomía a la tabla
sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(psd5)[rownames(sigtab), ], "matrix"))

# Manipulaciones varias para finalmente graficar los resultados
sigtabgen = subset(sigtab, !is.na(Genus))
# Phylum order
x = tapply(sigtabgen$log2FoldChange, sigtabgen$Phylum, function(x) max(x))
x = sort(x, TRUE)
sigtabgen$Phylum = factor(as.character(sigtabgen$Phylum), levels=names(x))
# Genus order
x = tapply(sigtabgen$log2FoldChange, sigtabgen$Genus, function(x) max(x))
x = sort(x, TRUE)
sigtabgen$Genus = factor(as.character(sigtabgen$Genus), levels=names(x))
ggplot(sigtabgen, aes(y=Genus, x=log2FoldChange, color=Phylum)) + 
    geom_vline(xintercept = 0.0, color = "gray", size = 0.5) +
    geom_point(size=4) + 
    theme(axis.text.x = element_text(angle = -90, hjust = 0, vjust=0.5, size = 10), axis.text.y = element_text(size = 13), legend.text = element_text(size = 13) )

df_res <- cbind(as(res.PR.sCR, "data.frame"), as(tax_table(psd5)[rownames(res.PR.sCR), ], "matrix"))
df_res <- df_res %>%
  rownames_to_column(var = "OTU") %>%
  arrange(padj)

(fdr_otu <- df_res %>%
  dplyr::filter(padj < 0.05))   
##       OTU  baseMean log2FoldChange    lfcSE      stat       pvalue         padj
## 1  ASV143  90.57919     -23.647542 2.798983 -8.448620 2.947668e-17 1.143695e-14
## 2   ASV78 138.27535     -24.242780 3.338224 -7.262179 3.809041e-13 3.986267e-11
## 3  ASV105 112.93822     -23.953855 3.303113 -7.251903 4.109554e-13 3.986267e-11
## 4  ASV110 108.77143     -23.903707 3.290335 -7.264825 3.735217e-13 3.986267e-11
## 5  ASV198  50.85653     -22.848421 3.389657 -6.740629 1.577022e-11 1.223769e-09
## 6  ASV334  31.80996     -22.194766 3.351546 -6.622247 3.537789e-11 2.287770e-09
## 7  ASV330  25.72342     -21.897927 3.388500 -6.462424 1.030386e-10 5.711284e-09
## 8   ASV48  21.06044     -21.618912 3.390043 -6.377180 1.803782e-10 7.063576e-09
## 9   ASV61  20.96650     -21.614137 3.390046 -6.375766 1.820509e-10 7.063576e-09
## 10 ASV276  21.57568     -21.653614 3.390027 -6.387446 1.686793e-10 7.063576e-09
## 11  ASV95  18.96324     -21.466058 3.390115 -6.331955 2.420742e-10 8.046069e-09
## 12 ASV173  27.80827     -21.431890 3.389883 -6.322310 2.576820e-10 8.046069e-09
## 13 ASV228  18.09508     -21.409920 3.390150 -6.315330 2.695848e-10 8.046069e-09
## 14 ASV299  14.98841     -21.146695 3.369857 -6.275251 3.490714e-10 9.674264e-09
## 15 ASV128  14.00293     -21.052789 3.390374 -6.209577 5.312756e-10 1.374233e-08
## 16 ASV165  26.28655     -20.722081 3.389912 -6.112868 9.785659e-10 2.373022e-08
## 17 ASV147  76.99358      -9.171426 2.113207 -4.340051 1.424494e-05 3.251197e-04
## 18 ASV171  65.75344      -8.943829 2.569298 -3.481040 4.994711e-04 1.076638e-02
## 19 ASV215  46.91935      -8.456888 2.448534 -3.453857 5.526295e-04 1.128528e-02
## 20 ASV259  35.65692      -8.061065 2.652187 -3.039403 2.370474e-03 4.598719e-02
##     Kingdom           Phylum          Class             Order
## 1  Bacteria Actinobacteriota Actinobacteria Bifidobacteriales
## 2  Bacteria Actinobacteriota Actinobacteria Bifidobacteriales
## 3  Bacteria Actinobacteriota Actinobacteria Bifidobacteriales
## 4  Bacteria Actinobacteriota Actinobacteria Bifidobacteriales
## 5  Bacteria     Bacteroidota    Bacteroidia     Bacteroidales
## 6  Bacteria       Firmicutes     Clostridia    Lachnospirales
## 7  Bacteria Actinobacteriota Coriobacteriia  Coriobacteriales
## 8  Bacteria       Firmicutes        Bacilli   Lactobacillales
## 9  Bacteria       Firmicutes        Bacilli   Lactobacillales
## 10 Bacteria       Firmicutes     Clostridia   Oscillospirales
## 11 Bacteria       Firmicutes        Bacilli   Lactobacillales
## 12 Bacteria     Bacteroidota    Bacteroidia     Bacteroidales
## 13 Bacteria     Bacteroidota    Bacteroidia     Bacteroidales
## 14 Bacteria       Firmicutes     Clostridia   Oscillospirales
## 15 Bacteria       Firmicutes        Bacilli   Lactobacillales
## 16 Bacteria     Bacteroidota    Bacteroidia     Bacteroidales
## 17 Bacteria       Firmicutes        Bacilli   Lactobacillales
## 18 Bacteria       Firmicutes        Bacilli   Lactobacillales
## 19 Bacteria       Firmicutes        Bacilli   Lactobacillales
## 20 Bacteria       Firmicutes     Clostridia    Lachnospirales
##                Family                        Genus
## 1  Bifidobacteriaceae              Bifidobacterium
## 2  Bifidobacteriaceae              Bifidobacterium
## 3  Bifidobacteriaceae              Bifidobacterium
## 4  Bifidobacteriaceae              Bifidobacterium
## 5      Bacteroidaceae                  Bacteroides
## 6     Lachnospiraceae                  Coprococcus
## 7     Eggerthellaceae             Senegalimassilia
## 8    Lactobacillaceae            Ligilactobacillus
## 9    Lactobacillaceae            Ligilactobacillus
## 10    Ruminococcaceae             Faecalibacterium
## 11   Lactobacillaceae            Ligilactobacillus
## 12     Bacteroidaceae                  Bacteroides
## 13     Bacteroidaceae                  Bacteroides
## 14    Ruminococcaceae             Faecalibacterium
## 15   Lactobacillaceae            Ligilactobacillus
## 16     Bacteroidaceae                  Bacteroides
## 17   Streptococcaceae                Streptococcus
## 18   Streptococcaceae                Streptococcus
## 19   Streptococcaceae                Streptococcus
## 20    Lachnospiraceae [Ruminococcus] torques group
#Metadata counts
as(sample_data(psd5), "data.frame") -> aaa
library(plyr)
counts <- ddply(aaa, .(aaa$Respuesta_100, aaa$Sexo), nrow)
names(counts) <- c("Respuesta_100", "Sexo", "Freq")
counts
##   Respuesta_100   Sexo Freq
## 1            PR HOMBRE    4
## 2            PR  MUJER    2
## 3           sCR HOMBRE   12
## 4           sCR  MUJER    5
list(fdr_otu$OTU)
## [[1]]
##  [1] "ASV143" "ASV78"  "ASV105" "ASV110" "ASV198" "ASV334" "ASV330" "ASV48" 
##  [9] "ASV61"  "ASV276" "ASV95"  "ASV173" "ASV228" "ASV299" "ASV128" "ASV165"
## [17] "ASV147" "ASV171" "ASV215" "ASV259"
#Create a phyloseq object OTUS significatives



GP.chl = prune_taxa(taxa_names(psd5) %in% c("ASV143", "ASV78",  "ASV105", "ASV110", "ASV198", "ASV334", "ASV330", "ASV48",  "ASV61","ASV276", "ASV95",  "ASV173", "ASV228", "ASV299", "ASV128", "ASV165", "ASV147", "ASV171", "ASV215", "ASV259"), psd5)

GP.chl
## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 20 taxa and 23 samples ]
## sample_data() Sample Data:       [ 23 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 20 taxa by 6 taxonomic ranks ]
## phy_tree()    Phylogenetic Tree: [ 20 tips and 19 internal nodes ]
taxa_sums(GP.chl)
##  ASV48  ASV61  ASV78  ASV95 ASV105 ASV110 ASV128 ASV143 ASV147 ASV165 ASV171 
##   3673   3338   2556   2609   2089   2027   1890   1698   1729   1885   1575 
## ASV173 ASV198 ASV215 ASV228 ASV259 ASV276 ASV299 ASV330 ASV334 
##   1849   1483   1072   1176    861    436    355    686    647
#Comparation of abundance
taxa_names(GP.chl)
##  [1] "ASV48"  "ASV61"  "ASV78"  "ASV95"  "ASV105" "ASV110" "ASV128" "ASV143"
##  [9] "ASV147" "ASV165" "ASV171" "ASV173" "ASV198" "ASV215" "ASV228" "ASV259"
## [17] "ASV276" "ASV299" "ASV330" "ASV334"
#metadata, taxonomíaa y abundancia del taxon más abundante de cada muestra

df2 <- psmelt(GP.chl)

head(df2) %>%
  kable(format = "html", col.names = colnames(df2)) %>%
  kable_styling() %>%
  kableExtra::scroll_box(width = "100%", height = "400px")
OTU Sample Abundance SampleID Edad_a_la_inclusion Sexo Respuesta_100 aloTPH_previo_si_no Citogenética_alto_riesgo Afect_Extramedular_PET_screening_si_no Dias_producción CRS_si_no Grado_máx_CRS HOSPITAL Visit Kingdom Phylum Class Order Family Genus
379 ASV48 HCB12-V2 1843 HCB12-V2 64 HOMBRE sCR NO SI NO 110 NO NA HCB V2 Bacteria Firmicutes Bacilli Lactobacillales Lactobacillaceae Ligilactobacillus
401 ASV61 HCB12-V2 1778 HCB12-V2 64 HOMBRE sCR NO SI NO 110 NO NA HCB V2 Bacteria Firmicutes Bacilli Lactobacillales Lactobacillaceae Ligilactobacillus
381 ASV48 HVA04-V2 1272 HVA04-V2 74 MUJER sCR NO NO NO 100 SI 10 HVA V2 Bacteria Firmicutes Bacilli Lactobacillales Lactobacillaceae Ligilactobacillus
442 ASV95 HCB12-V2 1115 HCB12-V2 64 HOMBRE sCR NO SI NO 110 NO NA HCB V2 Bacteria Firmicutes Bacilli Lactobacillales Lactobacillaceae Ligilactobacillus
134 ASV165 ISBS03-V2 1061 ISBS03-V2 70 MUJER sCR NO SI 100 SI 10 ISBS V2 Bacteria Bacteroidota Bacteroidia Bacteroidales Bacteroidaceae Bacteroides
399 ASV61 HVA04-V2 1029 HVA04-V2 74 MUJER sCR NO NO NO 100 SI 10 HVA V2 Bacteria Firmicutes Bacilli Lactobacillales Lactobacillaceae Ligilactobacillus
#Filtrar OTUS sifnificativas
xx <- df %>%
  group_by(OTU)%>%
  filter(OTU %in% c("ASV143", "ASV78",  "ASV105", "ASV110", "ASV198", "ASV334", "ASV330", "ASV48",  "ASV61","ASV276", "ASV95",  "ASV173", "ASV228", "ASV299", "ASV128", "ASV165", "ASV147", "ASV171", "ASV215", "ASV259"))

xx
## # A tibble: 460 × 21
## # Groups:   OTU [20]
##    OTU    Sample    Abundance SampleID  Edad_a_la_inclusion Sexo   Respuesta_100
##    <chr>  <chr>         <int> <chr>                   <int> <chr>  <fct>        
##  1 ASV48  HCB12-V2       1843 HCB12-V2                   64 HOMBRE sCR          
##  2 ASV61  HCB12-V2       1778 HCB12-V2                   64 HOMBRE sCR          
##  3 ASV48  HVA04-V2       1272 HVA04-V2                   74 MUJER  sCR          
##  4 ASV95  HCB12-V2       1115 HCB12-V2                   64 HOMBRE sCR          
##  5 ASV165 ISBS03-V2      1061 ISBS03-V2                  70 MUJER  sCR          
##  6 ASV61  HVA04-V2       1029 HVA04-V2                   74 MUJER  sCR          
##  7 ASV95  HVA04-V2       1020 HVA04-V2                   74 MUJER  sCR          
##  8 ASV173 ISBS03-V2      1003 ISBS03-V2                  70 MUJER  sCR          
##  9 ASV198 ISBS03-V2       819 ISBS03-V2                  70 MUJER  sCR          
## 10 ASV128 HCB12-V2        780 HCB12-V2                   64 HOMBRE sCR          
## # … with 450 more rows, and 14 more variables: aloTPH_previo_si_no <chr>,
## #   Citogenética_alto_riesgo <chr>,
## #   Afect_Extramedular_PET_screening_si_no <chr>, Dias_producción <int>,
## #   CRS_si_no <chr>, Grado_máx_CRS <int>, HOSPITAL <chr>, Visit <chr>,
## #   Kingdom <chr>, Phylum <chr>, Class <chr>, Order <chr>, Family <chr>,
## #   Genus <chr>
#Save table in .csv file
write_csv(xx, file="xx.csv")